Probing thermality beyond the diagonal

We investigate the off-diagonal sector of eigenstate thermalization using both local and non-local probes in 2-dimensional conformal field theories. A novel analysis of the asymptotics of OPE coefficients via the modular bootstrap is performed to extract the behaviour of the off-diagonal matrix elements. We also probe this sector using semi-classical heavy-light Virasoro blocks. The results demonstrate signatures of thermality and confirms the entropic suppression of the off-diagonal elements as necessitated by the eigenstate thermalization hypothesis.


Introduction
The interplay of ordering and fluctuations is important to understand the origin of statistical laws governing phases of matter. In the thermodynamic limit these laws begin to operate and control collective behaviour with an accuracy increasing with the number of degrees of freedom. Most scenarios involve a coarse-graining that gives rise to the exact statistical laws. The fluctuations are inherently quantum mechanical, and as we coarse grain over increasing number of quantum states, statistical thermodynamics takes over and just a few parameters like temperature, pressure, and volume describes the system. This is especially true if the system possesses a finite energy density that we can measure. When this happens we say that the system under consideration has thermalized. Thermal states in quantum mechanics are mixed states and there is no way that unitary operations can get us from a pure state to a thermal one. Naturally, this gives rise to an observable dependent notion of thermalization, i.e., one can look at the coarse grained expectation value of an operator, ψ ψ|O|ψ and try to approximate it with a thermal expectation value Z(β) −1 Tr[O e −βH ]. When the operator O is arbitrary and has support over a certain subregion V , then the coarse-grained reduced density matrix is well approximated by the reduced thermal density matrix. Usually the quantum states under consideration will have a complicated time evolution and consequently it is interesting to understand how it approaches thermalization at late times.
The eigenstate thermalization hypothesis (ETH) provides a mechanism (or, strictly speaking, a criterion) for closed quantum systems to be well described at late times by thermal equilibrium under unitary evolution [1][2][3]. The idea centers around the notion that Hamiltionian time evolution of eigenstates is trivial and, therefore, finite energy-density eigenstates are approximately thermal. The statement of the ETH proposes the following ansatz for the matrix elements of appropriate few-body observables (1.1) Here, E avg = (E m + E n )/2, S(E avg ) is the thermal entropy at the average energy E avg and ω = E m − E n . The functions g O and f O are smooth functions of their arguments. R mn denotes random variables distributed with a zero mean and unit variance. The connection of this ansatz with statistical mechanics is through g O (E) which is the statistical mechanics prediction for the operator O at the mean energy E. Eigenstate thermalization in this manner has been observed on the lattice for many non-integrable models, e.g. [3][4][5][6][7][8][9][10]. A crucial component of the above ansatz is the requirement that the off-diagonal elements are exponentially suppressed compared the diagonal ones. The intuition behind this proposed suppression can be illustrated as follows [11]. Consider a system in a pure state which is a generic superposition of the energy eigenstates, |ψ = m C m |m . The unitary evolution of an observable in this state is simply given by Let us make a rough estimate on these fluctuations. For an appropriately normalized pure state built from a superposition of N eigenstates, we have C m ∼ 1/ √ N . The fluctuation is therefore In the second step we have used the fact that the average amplitude of M random phases grows as √ M [12]. This shows that the typical values of the off-diagonal elements determine the fluctuations at late times. A similar estimate shows that the typical values of diagonal elements determine the time-average O(t) ∼ m|O|m typical . Therefore, in order for the fluctuations to be small to allow equilibriation, a necessary requirement is that the values of typical off-diagonal matrix elements should be much smaller than the diagonal ones. 1 This suppression is also suggested from close relations between quantum chaos and thermalization. The chaotic behaviour is effectively described by random matrix theory which also predicts the precise factor of e −S/2 for off-diagonal elements of hermitian operators. This suppression has been verified using numerical simulations in a number of lattice models [15] (including the SYK [16,17]). However, there has been a lack of analytical handle so far to confirm this prediction. This can be partly attributed to the lack of integrability in systems with a random matrix like behaviour. In this work, we take some steps in this direction by focussing on two dimensional conformal field theories. In 2d CFTs, the validity and related consequences which follow from (1.1) have been investigated in some works [18][19][20].
The constraints from modular invariance allow us to extract features of both the diagonal as well as the off-diagonal matrix elements. In the CFT parlance, the matrix elements of local observables are nothing but the OPE coefficients. The modular bootstrap that we shall perform uses the two-point function on the torus as the starting point. We shall show that the off-diagonal elements, when coarse-grained , are suppressed at least by the factor of e −S/2 , thereby providing evidence to the expectation (1.1) for a typical matrix element. This analysis is along the lines of the recent progress on extraction of asymptotics of OPE coefficients using modular features of correlation functions [21][22][23]. We shall also utilize statistics of the OPE coefficients and the inequality between the 1-norm and the 2-norm, to provide a new constraint on the light data of the CFT.
We also examine the off-diagonal sector using a bi-local probe of two light operators L(x)L(y). This observable has been intensely explored by a variety of approaches in the context of eigenstate thermalization, holography and black holes [24][25][26][27][28][29][30][31][32][33][34][35][36][37][38]. We shall study the off-diagonal elements of this probe at large central charge using the monodromy method for conformal blocks. For the natural choice of perturbation parameters, in the light operator dimensions, the exchanged operator dimension and in the difference in dimensions of the heavy primaries (H 1 and H 2 ), the off-diagonal conformal blocks display thermal features. The inverse temperature now is given by β = L(12∆ avg /c − 1) −1/2 , where ∆ avg = (∆ H 1 + ∆ H 2 )/2. This is the expected temperature consistent with the ETH ansatz, similar to the diagonal case [26]. This observable has also been analysed beyond the diagonal sector previously [28]. As we shall show, the blocks from the monodromy method cover a different regime in the parameter space. The agreement in the large c limit and deviations away from it will also be seen by comparing with blocks obtained via the Zamolodchikov recursion relations [39]. A dual holographic version of these conformal blocks will also be discussed.
The outline of this paper is as follows. In section 2 we extract the asymptotics of mean-squared OPE coefficients using modular bootstrap of torus 2-point functions. The asymptotics will be derived for generic states first and then refined for primaries using properties of torus blocks. In the following section 3 we discuss the statistics of our results and, in particular, use it together with previous results from [21] to obtain a constraint on OPE coefficients involving light excitations. Section 4 contains the analysis of off-diagonal elements of observables L(x)L(y) using the monodromy method and comparison to other approaches. We conclude in section 5.
2 Off-diagonal one point functions

Modular properties of correlation functions on the torus
Conformal field theories in 2-dimensions can be uniquely specified by their central charge c, their spectrum, i.e., a set of conformal weights {(h i ,h i )} of primaries O i , and coefficients C ijk appearing in 3-point correlation functions. The coefficients are also directly related to the fusion coefficients C k ij which tell how the field O k contributes to the operator product expansion of the primaries O i and O j . Any correlation function can be constructed in terms of these data. Bootstrap is a strategy to gain information about the above data by using various kinds of consistency relations. The most prominent bootstrap method comes from constraints imposed by crossing symmetry of four-point correlation functions on the sphere, which has led to impressive results in higher-dimensional CFTs [40][41][42][43][44][45]. Another type of constraint special to CFT in 2d comes from modular covariance. Defining CFTs consistently on the torus implies specific transformation properties of correlation functions under modular transformations. In particular, the partition function must be modular invariant, which gives highly non-trivial constraints on the spectrum. For example, these allow the ADE classification of minimal models [46] and also imposes universality of the high energy density of states in any (unitary) 2d CFT [47].

Set-up
Consider a torus specified by the modular parameter τ . All modular transformations, i.e., transformations that lead to equivalent tori, can be generated by the two transformations The group of modular transformations is SL(2, Z) and the most general transformation of the modular parameter is γ · τ = aτ +b cτ +d . We shall be interested in the correlation function of primary operators on the torus. The operator is located at the elliptic variable, w, which transforms to w/(cτ + d). Therefore, primaries transform as Two-point correlation functions on the torus are defined as with, q = e 2πiτ . This quantity is doubly-periodic in w 12 . Note that this is the unnormalized two-point function on the torus in which we do not divide by the partition function. Using the definition (2.3), equation (2.2) along with the fact that the partition function is modular invariant, it can be seen that the correlator has the following modular transformation The torus, T 2 ≡ S 1 β × S 1 L , describes a CFT on a spatial circle L and at finite temperature β. The modular parameter τ = iβ/L is purely imaginary. It follows immediately from (2.4) that there is a direct relation between the high and low temperature behaviour of the 2point functions via the S-modular transformation which takes τ → −1/τ . In 2d CFTs the low temperature behaviour is dictated by the light spectrum of the theory and allows very good approximations on various quantities. For example, in case of the partition function one can directly relate the asymptotic density of high energy states to the energy of the ground state (Casimir energy) which in turn is determined by the central charge [47]. In [21] a similar analysis was carried out using torus one point functions of primary operators. This yielded an average value of three point coefficients C EOE , where E labels an operator H with large dimensionless energy and O denotes a primary operator which is light compared to H. The average is over all states with that high energy. 2 We shall now apply the modular bootstrap method on the torus two-point functions. This will lead to an average value for C 2 EOE for two high energy operators H at energy E and H at energy E , and some scalar primary operator O.

Asymptotics of OPE coefficients
We begin with the thermal two-point function of a primary operator O. Without loss of generality, we place the operators at (0, t) and (0, 0) where H = 2π L (L 0 +L 0 − c 12 ). At low temperatures (L/β → 0) the leading contribution comes from the vacuum and this becomes the two-point function on the cylinder of circumference L which is completely fixed by conformal invariance Since we are considering the unnormalized correlator, we have the factor of Z(β) = e πcβ 6L . The subleading corrections are contributions from non-vacuum states which are exponentially suppressed as e −β∆χ . Here ∆ χ is the conformal dimension of the lightest operator in the CFT for which the four-point function χOOχ is non-zero -this may be a primary or the lightest descendant of the vacuum.
Using the S-modular transformation of the two point function on the torus we are led to the high temperature result (β/L → 0) (2.7) The terms in the ellipsis are suppressed in powers of e −2πL/β . Note that in general we consider complexified time in this text, where periodicity in the real part corresponds to the thermality of the system (i.e., the KMS condition) and evolution in physical time t corresponds to the imaginary part 3 Expanding the RHS of (2.5) with an insertion of a complete set of states then gives ( Here E i is the dimensionless energy (= conformal dimension ∆ i ) of the eigenstate |i which includes all states of the CFT, both primaries and their descendants. An integral representation of the above formula is where we have introduced the weighted spectral density (ω = E − E) (2.10) It is clear that (2.9) has the structure of a Laplace transform L in ∆ and a Fourier transform F in ω of the weighted spectral density, i.e., The weighted spectral density can then be obtained by inverting this expression. We get In the second line we have switched to the dimensionless variables,β = 2πβ/L andt = 2πt/L. The result for the Fourier transform in the third line follows from Mellin-Barnes integrals (see e.g., [49] and appendix A). We have taken ω > 0, therefore we closed the t-contour in the lower half plane, hence picked the −i prescription in (A.2) for the retarded 2-point function. Note that the combination E + ω/2 appearing in (2.12) is the average To perform the inverse Laplace transform we may utilize the saddle point approximation. Considering just the exponential factor in the integrand (2.12), we find that the saddle is located at .

(2.13)
This leads to an 'effective temperature' which is given in terms of the average (E + ω/2) of the dimensionless energies of the (off-)diagonal states. This thermal feature is exactly the same as observed earlier for the case of two point function of two light operators [28]. In [28] the same parameter was used to perform a uniformization transformation to calculate (off-)diagonal heavy light blocks. To justify the saddle point approximation, we will evaluate (2.12) by keeping all the factors in the integrand. In what follows we will keep ω fixed while taking E → ∞. The states with large E are expected to dominate at high temperatures, β 1. Hence setting x =βω, we can do an expansion around x → 0. (2.14) Using the integral representation of the modified Bessel function one can perform the integral (2.12) within the sum to get Since we are considering only large E asymptotics, we expand the Bessel function for large arguments. At large z, I ν (z) ≈ e z / √ 2πz, which in particular is independent of ν. This allows us to resum the series to recover the Γ functions, resulting finally in L 2π Here ρ(E) is the asymptotic density of high energy states given by Cardy's formula [47] ρ is the entropy devoid of any logarithmic corrections. The result (2.17) agrees with the saddle point result of (2.12). However, since we actually work with a contour integral in the complex plane, we prefer the method presented. We can now provide the asymptotics of the mean squared matrix element, |C EOE | 2 . It can be computed from where we can use (2.18) only for both E and E large. The average here is over all heavy states of the CFT and, we reiterate that it does not distinguish between primaries and descendants. Therefore, for the leading order approximation, the mean square OPE coefficient can be written as In the final step we have worked in the approximation |ω| = |E − E | is much smaller than E avg and also omitted all sub-leading contributions for large E avg . This approximation is justified since the Γ-function factor is peaked around ω = 0. The factor (12E avg /c − 1) ∆ O +1/4 above can be regarded as a logarithmic correction to the entropy in presence of the probe operators. The result demonstrates that mean square OPE coefficients of the above kind are entropically suppressed as e −S(Eavg) . The above equation (2.20) is one of the main results of this work.
It is worthwhile to interpret the result (2.20) holographically. The quantity |C EOE | 2 provides a measure of the transition rate of a black hole microstate of energy E to another of energy E along with the emission of a scalar (dual to the primary O). The factor of e −S is the probability of choosing a single black hole microstate and the Γ-function factor is the black hole emission rate for scalars [50] 4 .

Fluctuations and eigenstate thermalization
A measure of fluctuations of the operator O in the eigenstate |n can be obtained from the following quantity; see [52, §6] for further details.
The final relation is true only for states which obey the ETH ansatz (1.1). Moreover, the function f O (E avg , ω) is peaked around ω = 0 and thereby the dominant contribution to sum in (2.21) is from the small window of states within the regime E m − E n 1.
The form of the summand (2.21) is clearly consistent with (2.20) obtained using modular bootstrap. Concretely, for |n being a typical high energy eigenstate we can plug in the asymptotic mean square average (2.20) for | n|O|m | 2 in the first line of (2.21). We obtain the final equality which is also predicted by ETH, i.e., the final relation of (2.21). Notably, the entropic suppression e −S(Em+ω/2) comes out precisely as expected. The average of the factor |R mn | 2 is a constant. Furthermore, we can read off the function f O (E avg , ω) This is an explicit verification of smoothness of the function f O as conjectured by the ETH ansatz (1.1). It also agrees with the prediction that this is a function in the average energy of the in and out states E avg and the frequency ω characterizing the off-diagonal elements. It is also crucial to observe that the distribution of the matrix elements with respect to ω is not that of a generic Gaussian ensemble. The distribution is actually that of generalized hyperbolic secants [53]. This should be contrasted with the cases where the underlying behaviour is believed to be governed by random matrices [3,15,54,55].

Refining asymptotics for primaries
The formula (2.20) is for the averaged matrix elements of high energy eigenstates regardless of them being primaries or descendants. In this subsection, we show that it can be refined The thermal Euclidean two point function of light primaries O located at elliptic coordinates (w,w) and (0, 0) on the torus T 2 admits the following expansion, involving a sum only over the exchanged primaries Here, q = e 2πiτ , z = e 2πiw and h ij = h i − h j . The factor G O ij (z, q) and its antiholomorphic counterpart are the 2-point torus blocks which encode the contributions of the descendants of primaries. Note that we have considered the expansion along the projection channel as illustrated in Fig. 1. This channel is also referred to the s-channel or the necklace channel literature. If we had chosen to work with the other channel (called the t-channel or the OPE channel) coefficients of the kind |C OO∆ i | 2 would have appeared in the sum above. As in the previous subsections, we shall work with the rectangular torus such that q = e −β (setting L = 2π for convenience) and the locations of the operators are at (0, t E ) and (0, 0). The 2-point torus blocks in the projection channel have the following form [56,57] and similarly for the anti-holomorphic part. B P |Q i,j appearing above are the elements of the inverse Gram matrices of the exchanged primaries ∆ i,j . All the inner products in the torus block (2.24) are completely determined by the Virasoro Ward identities.
We are interested in the high temperature regime of the two-point function (2.23) and we specialize to the situation where the operators O are light. In this regime, the dominant contribution arises from heavy primaries. We shall also assume h ij h i . The torus blocks for the heavy exchanged primaries have a greatly simplified form. In order to see this, note that the elements of the inverse Gram matrices in the regime h i , h j 1 always have a leading behaviour as some inverse power of h i , h j . On the other hand, the inner products in the numerator, ∆ i , P |O(z)|Q, ∆ j , behave as some positive power of h i for h i h O 5 and h i h ij . It can be seen that only if M = S and T = N in both inner products in (2.24), we get an order one contribution from the summand. Otherwise, these get suppressed by inverse powers of h i , h j . The net effect due to the presence of the inverse Gram matrices is, therefore, to restrict the dominant contribution to the 'diagonal sector' only The sum over N is over all possible descendants at each level n which can be uniquely labelled by partitions of integers. Retaining just the leading terms (which are order-one) in the h i , h j → ∞ limit, merely counts the number of descendants at each level This is nothing but the character of the primary h i (together with the prefactor accounted for in (2.23)). To summarize in words, the heavy-light regime of the torus blocks can be approximated by the character of the heavy primary, i.e., the deformations caused by the presence of the light operators are negligible. This is analogous to the situation for the lower-point torus block in [21,58,59]. This can be checked to a very low order in the q-expansion using the results of [57, Appendix A].
Using the S-modular transformation of the Dedekind-eta function, η(− 1 τ ) = √ −iτ η(τ ), we arrive at the high temperature (β → 0) behaviour of the following product The corrections to the above are in inverse powers of h i,j and additional exponentially suppressed corrections inq = e −4π 2 /β . We can now in a position to define the weighted spectral density over primaries As before, the high temperature limit of the Euclidean thermal two point function is 3β . (2.29) Using (2.27) and (2.29) and the high-energy/high-temperature version of (2.23) we have Next we analytically continue both sides to Lorentzian time and follow the same steps as in the previous section to solve for J P O by performing an inverse Fourier transform followed by the inverse Laplace transform. We keep only the leading order contributions for high energies. Our final result for the weighted spectral density of heavy primaries then is (with We have used the density of heavy primaries [21] ρ Using the latter we can also write where |C ∆O∆ | 2 now is the mean squared OPE coefficients with the average taken over all heavy primaries with conformal dimension ∆ and ∆ . The final result of this section therefore is where again the last line is true when |ω| is much smaller than ∆ avg and we define The form of the above average is similar to that of (2.20). The precise shift in the central charge from c to c − 1 has also been observed in mean heavyheavy-light OPE coefficient in [21]. This is possibly related to the fact that the torus block (2.26) in the heavy limit (or the characters of primaries in c > 1 Virasoro CFTs, which count the descendants) is the same as the partition function of the c = 1 free boson theory. It would be interesting to derive the above result holographically building on some of the techniques of [57,59,60].

A bound from the statistics of OPE coefficients
The results of the previous section using modular bootstrap provide statistics of high energy CFT data. The expectation values of a randomly chosen probe are largely dependent on the actual distribution of the quantity of interest. A good measure for this is the squared variance σ 2 = X 2 − X 2 , which quantify how peaked a distribution on a random variable X is around its expectation value X . The squared variance is always positive, s.t. it follows immediately that our result on the mean squared OPE coefficient sets an upper bound for the mean OPE coefficients of [21] E|O|E Note that the average here is over two sets of heavy eigenstates which are specified by E and E . Since |C EOE | 2 is entropically suppressed, so is the variance. The distribution of OPE coefficients C EOE is therefore rather sharply peaked around its mean value. This justifies the notion of a typical (∼ almost all) OPE coefficient for which one gets This shows that in a 2d CFT typical off-diagonal OPE coefficients of high energy eigenstates are at least suppressed as expected from ETH ansatz (1.1). Our considerations here, however, do not allow provide a lower bound.
If we consider two different heavy states but with the same energy, i.e. ω = 0, the mean squared OPE coefficient takes the form Note that the average here includes off-diagonal OPE coefficients with different states but having the same energy. This situation is possible since the heavy spectrum is highly degenerate. To make this clear we shall slightly change the notation. The average coming from purely diagonal entries (i.e., when both the energy and the state are the same) is bounded from above as 6 This upper bound is in particular not entropically suppressed. In the diagonal case (energy and state) the mean OPE coefficient was computed in [21] from the modular properties of one point functions on the torus. The result is where χ is the lightest field in the spectrum of the CFT with C χOχ = 0. Note that also in this case (3.1) must hold which together with (3.5) is consistent with all previous assumptions made in this text and in [21], in particular with E 1 and the unitarity bound, ∆ χ > 0. We can even go further and constrain C χOχ when ∆ χ is very small. Using (3.1) with (3.5) we arrive at the following inequality (3.7) On the RHS the central charge c and the operator dimension ∆ O are fixed by choice, and ∆ χ is fixed by the assumption that χ is the lightest field in the spectrum of the CFT that has non-vanishing OPE coefficient C χOχ . The inequality is true for any large enough E. In fact, the RHS has a minimum at which is large for ∆ 2 χ c 192π 2 . This implies that if the latter condition on the dimension of χ is met, we can constrain the OPE coefficient C χOχ (by plugging in E = E min in (3.7)) This provides a novel universal constraint on light data of a CFT. As it can be seen, the derivation uses modular bootstrap and the statistics that follow from it. This is a surprising by-product of our analysis.

Off-diagonal two-point functions
The analysis of the previous section provides results of the averaged OPE coefficients which are consistent with the ETH ansatz. In this section we shall work with a different probe and investigate the off-diagonal elements thereof. In what follows, the results fall under the context of the stronger form of ETH which posits that single energy eigenstates are thermal [61]. The following calculation will however be confined to the large central charge regime.
The perturbative monodromy method can be used to derive an off-diagonal generalization of the heavy-light conformal block V P (h H 1 , h H 2 , h L , c; z) in the large c limit with all the scaled operator dimensions h i /c held fixed and set hierarchically. The heavy conformal dimensions (h H 1 , h H 2 ) much larger than the light ones (h L ) and the exchanged primary (h P ). We note that, in a different regime of parameter space, with h H 1 /c, h H 2 /c, h H 1 −h H 2 , h P , h L fixed and c → ∞, V P has been obtained from the global conformal block using a background field method in [28].

Off-diagonal conformal blocks using the monodromy method
We consider the correlation function two heavy and two light operators. The conformal partial wave expansion of this correlator is, (4.1) Here, H 1 and H 2 are heavy operators and, in general, H 1 = H 2 . Here we choose H 1 > H 2 without loss of generality. Note that by conformal transformations we could have put the insertions at (x 1 , x 2 , x 3 , x 4 ) = (0, x, 1, ∞). However, since we shall like to make contact with a holographic interpretation it turns out to be transparent to keep the light operators at x and y. We shall fix our attention to the holomorphic part of the conformal block, V P (h i , c; x, y). In the asymptotic heavy limit the states correspond to black hole geometry in the bulk where the back-reaction due to light operators is negligible. As mentioned earlier, the monodromy technique that we use [24] 7 is valid only in the large central charge limit, whilst keeping the ratios i ≡ 6h i /c fixed. In this limit, the blocks exponentiate The monodromy method considers a conformal block of interest with the presence of an auxiliary fieldψ(z) which has a null state at level 2. The null state condition can be translated into the differential equation for ψ(z, x, y) ≡ H 1 (z 1 )L(z 2 )L(z 3 )|p p|ψ(z)H 2 (z 4 ) . Here T (z, z i ) is the stress-tensor wavefunction in the presence of the heavy and light operator insertions. The strategy is to solve for the function f P in (4.2) by demanding consistent monodromies of the solutions of the above differential equation order by order in perturbation theory. The perturbative parameters in our case is in L , P and H 1 − H 2 ≡ 2 . Thus in addition to the exponentiation of block limits, we also work in the limit when, ( H 1 + H 2 )/2 ≡¯ L , , P . The c i are accessory parameters, c i = ∂ z i f P . We also have the behaviour T (z) ∼ z −4 for large z, which imposes the constraints These constraints allow us to re-express the stress-tensor wavefunction as, . (4.5) As mentioned earlier, in the heavy-light limit¯ is larger than other parameters, hence the zeroth order homogenous equation (4.3) is which has two solutions ψ The zeroth order solution can be used to find the perturbative solution at the first order, ψ (1) ± , using the method of variation of parameters. This solution comes with a non-trivial monodromy. If we encricle ψ This equation can be solved for c x which is given by where we defined, u = y/x. Note that in this case we have two non-trivial accessory parameters, c x = ∂ x f P and c y = ∂ y f P . The constraints (4.4) yield y∂ y f P + x∂ x f P = 2( + L ). (4.9) The general solution of the above is Now differentiating the above with respect to log x we obtain x∂ x f P = 2( + L ) − ug (u). (4.11) This object by definition is just xc x which is a function of u (4.8). Therefore we can now integrate (4.11) to find g(u). The root in (4.8) and the constant of integration is fixed so as to yield the expected OPE limit of the conformal block. Putting together everything, we can bring the result to the following form (till the linear order , L , P ) f H 1 ,H 2 ,L P (x, y) = log xy + 2 L log x αavg − y αavg α avg (xy) αavg−1 2 + P log α avg (x αavg/2 + y αavg/2 ) 4(x αavg/2 − y αavg/2 ) . (4.12) In the diagonal limit, = 0 and the operator locations x → 1 and y → x the answer reduces to [32, equation (2.26)]. The conformal block can then be obtained from the above function by exponentiation (4.2) which is the main result of this section. For a sparse spectrum of light operators, this result provides an analytic expression of the smooth function f LL (¯ , ) which appears in the version of the ETH ansatz (1.1) for the bilocal probe L(x)L(y). We notice that there is a power law decay (xy) −(h H 1 −h H 2 )/2 for the off-diagonals. The full 4-point function also appears with the OPE coefficients C LLχ C H 1 H 2 χ . Here, χ is the lightest primary appearing in the fusion of H 1 and H 2 and two L operators. For a typical state, our estimate of the C H 1 H 2 χ from Section 2 shows the presence of entropic suppression also in this case. Note that, unlike the diagonal heavy-light correlator (in a typical holographic CFT), the dominant contribution is not from the vacuum block in the off-diagonal case. If the light spectrum is sufficiently sparse, the dominant contribution to H 2 LLH 1 is from the block of the light primary χ.
The parameter α avg becomes imaginary when (h H 1 + h H 2 )/2 > c/24. In the bulk dual this corresponds to the black hole regime. The conformal block (4.13) then acquires periodicities akin to thermal correlator. In particular, the second factor [· · · ] 2h L in (4.13) has the form of [sinh(π /β)] −2h L , when written in cylinder coordinates. Here, β is exactly the same as the saddle which we found in the modular bootstrap analysis -equation (2.13).

Zamolodchikov recursion
The Virasoro blocks can be computed at any value of central charge but are not known in closed form [39,62]. Nevertheless, the coefficients of the block in a small cross-ratio expansion can be determined algorithmically using Zamolodichikov's recursion 8 .
The S-channel (or z → 1) expansion of log V hp (c, h L , h H 1 , h H 2 , z) can be organized in the following manner In the regime of parameter space of the monodromy analysis, we have the following expansions for the coefficients v n (4.15) Note that the above expansion is obtained by performing a 1/c expansion first. This is followed by a scaling and then an expansion to the linear order in δ. These coefficients match exactly with those of the monodromy method (4.13).
It is worthwhile to observe that there are divergences arising from thermal periodicities in the conformal block calculated using the monodromy method. This happens for h H 1 + h H 2 > c/12 for which α avg is imaginary. These divergences are often referred to as forbidden singularities. As shown in [63], this is purely an artefact of the large central charge limit 9 . The blocks computed numerically using the Zamolodchikov recursion relations do not display these divergences. We have checked that the story is the same for the off-diagonal blocks.

Heavy-light blocks from the background field method
At large central charge, the heavy light conformal block has also been investigated in the regime in which the dimensions of the light and intermediate operators and differences between dimensions of heavy operators held fixed [28]. This is a different regime in the parameter space of conformal dimensions from the one covered by the monodromy method. The monodromy method has the ratios of the conformal dimensions and the central charge held fixed. The result of [28] is which also probes the off-diagonal sector. It can be shown to agree with the expansion using the Zamolodchikov recursion but using a different sequence of limits. One needs to expand first in h L , h p and h H 1 − h H 2 and then expand in 1/c. This result also remarkably agrees on the nose with a bulk computation [31].
It is interesting to contrast the block (4.16) with the monodromy block (4.13). In order to do this, we use the following Euler identity 10 for the hypergeometric factor Plugging this in (4.16) this gives The factors on the first line precisely agree with the h p -independent piece of the monodromy block (4.13). The piece that is dependent on the intermediate conformal dimension differs. Both the results agree for diagonal and vacuum case, h H 1 = h H 2 and h p = 0.

Holography
The result for the conformal block using the monodromy method (4.12) can be interpreted in the bulk dual via a configuration of geodesics on a fixed time slice of the global AdS 3 in presence of a conical singularity at the origin [32]. The metric is given by  and h H 2 . The dual CFT then lives on a cylinder.
The operators dual to the light scalar L and the one that labels the block P are shown in Fig. 2 by the green curve and the blue dotted line respectively. Both of them correspond to massive scalar fields, φ L and φ P in the bulk with mass given by m L = 2 h L (h L − 1) and similarly for m P . The last term in equation (4.12) is precisely the geodesic length of φ P , while the second term is the geodesic length of φ L times conformal factors coming from the transformation from the cylinder to the plane. The remaining term (i.e., the first one in (4.12)) is the product of regularized radial geodesic lengths, from the boundary anchoring points of φ L to the centre. 11 In the block both these come weighted with the same power of h H 1 − h H 2 . This is interpreted as yet another bulk massive scalar field, φ which is dual to an operator carrying conformal dimension given by the difference of the dimensions of the heavy operators. 12

Conclusions
In this work we have estimated the asymptotic behaviour of averaged matrix elements of an arbitrary scalar operator in the energy eigenbasis in the context of 2d CFTs. The offdiagonal matrix elements are universal and consistent with the eigenstate thermalization hypothesis. They dictate quantum fluctuations and, hence, play a key role in the investigation of non-equilibrium physics. Our analysis, which is valid for any positive value of the central charge, involves a coarse graining over exponentially large number of states. The final result for | E|O|E | 2 for large E, E with E −E small and non-zero is of the form e −S Γ, where S is the entropy at the average energy and Γ, which is written in terms of gamma functions, measures the spread of the fluctuations. We have also found the asymptotics of actual OPE coefficients corresponding to primaries using properties of 2-point torus blocks. These findings from modular bootstrap also provide information that allows to investigate late time behaviour of thermal 2-point functions on the torus. This can potentially furnish a 2d CFT analogue of a part of the analysis performed for the SYK model in [55, §7].
We then used our results along with results from [21] to provide a constraint using positivity of the variance of diagonal matrix elements. While we were able to give an upper bound to the variance of the diagonal matrix elements, it may be possible to find the diagonal variance itself. This will require estimation of a weighted sum of the form where the prime indicates that we sum over nonidentical states which are degenerate in their energies.
In the latter half of this work we studied off-diagonal Virasoro blocks of the form, H 1 LLH 2 in the large central charge limit. The calculation was done using the monodromy method where the conformal dimensions were taken to be of the order of the central charge. We performed the calculation with the conformal dimensions of both H 1 and H 2 much larger than those of P (intermediate exchange operator), L, and the difference in the dimensions of H 1 and H 2 . We obtained the answer for the block till linear order in h L /c, (h H 1 − h H 2 )/c and h P /c. The answer for the blocks bears some thermal features. Correspondingly, the holographic interpretation involves geodesics in a conical defect geometry. The scalar probes forming the geodesic network have masses proportional to h L , h P and h H 1 − h H 2 . The conical defect geometry can be analytically continued to a BTZ black hole geometry whose temperature is once again dictated by the average, (h H 1 + h H 2 )/2. In this case the thermal description emerges, not by coarse-graining, but at large central charge, i.e., in semi-classical limit.
There are various possible generalizations and applications of this work. It is interesting to look at the off-diagonal out of time ordered correlator of the form, W (t)V W (t)V β in a thermal state, when the conformal dimensions of W and W are only slightly different. Upon following the relevant analytic continuation prescription of the Euclidean correlator [64] one obtains for the OTOC, The dual interpretation is that of correlators of the light scalars dual to V in a shockwave background corresponding to the averaged weight of W and W . It would be interesting to investigate further implications of this. Furthermore, it would be worthwhile to explore a generalization of black hole collapse with a non-uniform distribution of matter by a continuum version of our monodromy calculation [65].
If ETH holds approximately in a local theory for all operators within a subregion A in the excited energy eigenstate, |ψ , then the reduced density matrix satisfies, ρ ψ A = Tr A c |ψ ψ| ≈ Tr A c e −βH e −βH A , where β is determined from H β = ψ|H|ψ and H A is the Hamiltonian restricted to the subregion A. The above equality can be used to approximately construct the full Hamiltonian given an eigenstate that satisfies ETH [61,66]. Once the conditions for validity of ETH for CFTs are well understood, the construction of a chaotic local Hamiltonian can be potentially possible, with a view towards building holographic CFTs from bottom-up.
If the conformal field theory has additional conserved currents, then the ensemble approximating of charged excited states is the grand canonical ensemble. The modular properties get modified when there is an additional Kac-Moody along with the Virasoro as the chiral algebra [67]. This gets reflected in the modular bootstrap analysis of diagonal OPE coefficients [22] by the appearance of spectral-flow invariants and Ahranov-Bohm like phases. It will be interesting to see how the presence of additional global symmetries modify the off-diagonal matrix elements and, in particular, how the fluctuations manifestly depend on the conserved charges. The work of SD is supported by the NCCR SwissMAP, funded by the Swiss National Science Foundation. SD also thanks AEI Potsdam and Swansea University for hospitality where parts of this work were completed. Finally, the authors thank the participants and organizers of 'Strings, Geometry & Black Holes' at King's College London for simulating discussions and an opportunity to present this work.

A Mellin-Barnes integrals
The following integral is useful for some of the calculations of this work. For any complex A and B The ±i prescription is chosen in order to avoid the singularities along the real t axis.