Conformal fields and operator product expansion in critical quantum spin chains

TL;DR

This paper introduces a variational method to identify lattice operators with conformal field theory operators in critical quantum spin chains, enabling numerical estimation of OPE coefficients, demonstrated on the critical Ising model.

Contribution

A novel variational approach to map lattice operators to CFT operators and compute OPE coefficients in critical quantum spin chains.

Findings

Successfully identified lattice primary operators in the critical Ising chain.

Numerically estimated the OPE coefficients matching CFT predictions.

Provided a practical method for connecting lattice models with conformal field theories.

Abstract

We propose a variational method for identifying lattice operators in a critical quantum spin chain with scaling operators in the underlying conformal field theory (CFT). In particular, this allows us to build a lattice version of the primary operators of the CFT, from which we can numerically estimate the operator product expansion coefficients $C_{\alpha\beta\gamma}^{\textrm{ CFT}}$. We demonstrate the approach with the critical Ising quantum spin chain.

Figures

Tables

Table 1. Table 1: Expansion ( 11 ) for simple lattice operators in the Ising model. Two-spin operators A ​ ( j ) ​ B ​ ( j + 1 ) 𝐴 𝑗 𝐵 𝑗 1 A(j)B(j+1) are organized into terms that are even or odd terms under exchange j ↔ j + 1 ↔ 𝑗 𝑗 1 j\leftrightarrow j+1 , e.g. X ​ Y ± Y ​ X plus-or-minus 𝑋 𝑌 𝑌 𝑋 XY\pm YX . The set 𝒜 𝒜 \mathcal{A} of CFT operators is given in Eqs. ( 17 )-( 19 ). Coefficients smaller than 5 × 10 − 3 5 superscript 10 3 5\times 10^{-3} are not shown. The number of significant digits is determined case by case by requiring that a particular digit does not change under extrapolation with different sets of system sizes up to N = 96 𝑁 96 N=96 Sup . Note that we omit the superscript CFT on the CFT scaling operators.

Lattice	CFT	Lattice	CFT
$X$	$0.803121\sigma -0.017{\partial }_{\tau }^2\sigma -0.033{\partial }_x^2\sigma$	$XY+YX$	$0.3184i{\partial }_{\tau }ϵ$
$Y$	$0.8031i{\partial }_{\tau }\sigma$	$XY-YX$	$0.6366\left(T-\overline{T}\right)$
$Z$	$0.636620𝟙-0.15915\left(T+\overline{T}\right)-0.31831ϵ+0.010{\partial }_{\tau }^2ϵ$	$XZ+ZX$	$0.803121\sigma -0.820{\partial }_{\tau }^2\sigma -0.736{\partial }_x^2\sigma$
$XX$	$0.636620𝟙-0.15915\left(T+\overline{T}\right)+0.31831ϵ-0.010{\partial }_{\tau }^2ϵ$	$XZ-ZX$	$1.205{\partial }_x\sigma$
$YY$	$-0.212207𝟙+0.4774\left(T+\overline{T}\right)+0.31831ϵ-0.089{\partial }_{\tau }^2ϵ$	$YZ+ZY$	$2.41i{\partial }_{\tau }\sigma$
$ZZ$	$0.540380𝟙-0.5403\left(T+\overline{T}\right)-0.54038ϵ+0.067{\partial }_{\tau }^2ϵ-0.051{\partial }_x^2ϵ$	$YZ-ZY$	$-0.4015i{\partial }_{\tau }{\partial }_x\sigma$
$XIX$	$0.540380𝟙-0.5403\left(T+\overline{T}\right)+0.54038ϵ-0.067{\partial }_{\tau }^2ϵ+0.051{\partial }_x^2ϵ$	$XZX$	$0.212207𝟙-0.4774\left(T+\overline{T}\right)+0.31831ϵ-0.089{\partial }_{\tau }^2ϵ$

Table 2. Table 2: Lattice operators with even Z 2 subscript 𝑍 2 Z_{2} symmetry and their representation using Majorana fermion operators

spin operator	fermion operator
$X\left(j\right)X\left(j+1\right)$	$-2i\psi \left(2j\right)\psi \left(2j+1\right)$
$Z\left(j\right)$	$-2i\psi \left(2j-1\right)\psi \left(2j\right)$
$X\left(j\right)Y\left(j+1\right)$	$-2i\psi \left(2j\right)\psi \left(2j+2\right)$
$Y\left(j\right)X\left(j+1\right)$	$2i\psi \left(2j-1\right)\psi \left(2j+1\right)$
$Y\left(j\right)Y\left(j+1\right)$	$2i\psi \left(2j-1\right)\psi \left(2j+2\right)$
$Z\left(j\right)Z\left(j+1\right)$	$-4\psi \left(2j-1\right)\psi \left(2j\right)\psi \left(2j+1\right)\psi \left(2j+2\right)$
$X\left(j\right)I\left(j+1\right)X\left(j+2\right)$	$-4\psi \left(2j\right)\psi \left(2j+1\right)\psi \left(2j+2\right)\psi \left(2j+3\right)$
$X\left(j\right)Z\left(j+1\right)X\left(j+2\right)$	$-2i\psi \left(2j\right)\psi \left(2j+3\right)$

Table 3. Table 3: Correspondence between lattice operators and CFT operators for the Ising model. The truncated set of CFT operators contains 𝟙 , T , T ¯ , ϵ , ∂ τ ϵ , ∂ τ 2 ϵ double-struck-𝟙 𝑇 ¯ 𝑇 italic-ϵ subscript 𝜏 italic-ϵ subscript superscript 2 𝜏 italic-ϵ \mathbb{1},T,\bar{T},\epsilon,\partial_{\tau}\epsilon,\partial^{2}_{\tau}\epsilon . Coefficients are obtained analytically in the top table. The bottom table is the same as the top table except that coefficients are shown their approximate values to 5 digits to compare with numerical results. The subscript ”CFT” is omitted in the column of CFT operators.

Lattice operator	CFT operator
$Z$	$\left(2/\pi \right)𝟙-1/\left(2\pi \right)\left(T+\overline{T}\right)-\left(1/\pi \right)ϵ+1/\left(32\pi \right){\partial }_{\tau }^2ϵ$
$XX$	$\left(2/\pi \right)𝟙-1/\left(2\pi \right)\left(T+\overline{T}\right)+\left(1/\pi \right)ϵ-1/\left(32\pi \right){\partial }_{\tau }^2ϵ$
$YY$	$-2/\left(3\pi \right)𝟙+3/\left(2\pi \right)\left(T+\overline{T}\right)+\left(1/\pi \right)ϵ-9/\left(32\pi \right){\partial }_{\tau }^2ϵ$
$ZZ$	$16/\left(3{\pi }^2\right)𝟙-16/\left(3{\pi }^2\right)\left(T+\overline{T}\right)-16/\left(3{\pi }^2\right)ϵ+2/\left(3{\pi }^2\right){\partial }_{\tau }^2ϵ$
$XZX$	$2/\left(3\pi \right)𝟙-3/\left(2\pi \right)\left(T+\overline{T}\right)+\left(1/\pi \right)ϵ-9/\left(32\pi \right){\partial }_{\tau }^2ϵ$
$XIX$	$16/\left(3{\pi }^2\right)𝟙-16/\left(3{\pi }^2\right)\left(T+\overline{T}\right)+16/\left(3{\pi }^2\right)ϵ-2/\left(3{\pi }^2\right){\partial }_{\tau }^2ϵ$
$-i\left(XY+YX\right)$	$\left(1/\pi \right){\partial }_{\tau }ϵ$
$XY-YX$	$\left(2/\pi \right)\left(T-\overline{T}\right)$

Equations

$$H^{\text{textsl}CFT}|{\psi }_{\alpha }^{\text{textsl}CFT}⟩$$

$$P^{\text{textsl}CFT}|{\psi }_{\alpha }^{\text{textsl}CFT}⟩$$

$$\mathcal{O}\sim \sum _{\alpha }{a}_{\alpha }{\psi }_{\alpha }^{\text{textsl}CFT}\equiv {\mathcal{O}}^{\text{textsl}CFT}.$$

$$C_{\alpha \beta \gamma }\equiv {\left(\frac{2\pi }{N}\right)}^{-{\Delta }_{\alpha }}⟨{\varphi }_{\beta }|{\varphi }_{\alpha }(0)|{\varphi }_{\gamma }⟩$$

$${\mathcal{O}}^{\text{textsl}CFT}(x)$$

$${\mathcal{O}}^{\text{textsl}CFT,s}$$

$${\varphi }_{\alpha }^{\text{textsl}CFT,s}|0^{\text{textsl}CFT}⟩={\left(\frac{2\pi }{L}\right)}^{{\Delta }_{\alpha }^{\text{textsl}CFT}}\sum _{m,\bar{m}\ge 0}{c}_{\alpha ,(m,\bar{m})}^s|{\varphi }_{\alpha ,(m,\bar{m})}^{\text{textsl}CFT}⟩,$$

$$c_{\alpha ,(m,\bar{m})}^s\equiv {\delta }_{m-\bar{m}+s_{\alpha }^{\text{textsl}CFT},s}\sqrt{\frac{\Gamma (m+2h_{\alpha }^{\text{textsl}CFT})\Gamma (\bar{m}+2{\bar{h}}_{\alpha }^{\text{textsl}CFT})}{m!\bar{m}!\Gamma (2h_{\alpha }^{\text{textsl}CFT})\Gamma (2{\bar{h}}_{\alpha }^{\text{textsl}CFT})}}.$$

$${\varphi }_{\alpha ,(k,\bar{k})}^{\text{textsl}CFT,s}|0^{\text{textsl}CFT}⟩={\left(\frac{2\pi }{L}\right)}^{{\Delta }_{\alpha }^{\text{textsl}CFT}+k+\bar{k}}\times$$

$$\displaystyle\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \sum_{m,\bar{m}\geq 0}(m+h^{\textsl{\tiny CFT}}_{\alpha})^{k}(\bar{m}+\bar{h}^{\textsl{\tiny CFT}}_{\alpha})^{\bar{k}}c^{s}_{\alpha,(m,\bar{m})}|\phi^{\textsl{\tiny CFT}}_{\alpha,(m,\bar{m})}\rangle.\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\$$

$$C_{\alpha \beta \gamma }^{\text{textsl}CFT}={\left(\frac{2\pi }{L}\right)}^{-{\Delta }_{\alpha }}⟨{\varphi }_{\gamma }^{\text{textsl}CFT}|{\varphi }_{\alpha }^{\text{textsl}CFT}(0)|{\varphi }_{\beta }^{\text{textsl}CFT}⟩,$$

$$\mathcal{O}\text{leavevmode}\stackrel{\text{mbox}approx}{\sim }\text{leavevmode}\sum _{\alpha \in \mathcal{A}}{a}_{\alpha }{\psi }_{\alpha }^{\text{textsl}CFT}\equiv {\tilde{\mathcal{O}}}^{\text{textsl}CFT}$$

$$⟨{\psi }_{\beta }|\mathcal{O}|{\psi }_{\alpha }⟩=⟨{\psi }_{\beta }^{\text{textsl}CFT}|{\tilde{\mathcal{O}}}^{\text{textsl}CFT}|{\psi }_{\alpha }^{\text{textsl}CFT}⟩+O\left(\frac{1}{N^{{\Delta }_c}}\right),$$

$$\mathcal{O}(j)=\sum _s{\mathcal{O}}^s{e}^{-is2\pi j/N},\text{leavevmode}\text{leavevmode}\text{leavevmode}\text{leavevmode}{\mathcal{O}}^s\equiv \frac{1}{N}\sum _{j=1}^N\mathcal{O}(j)e^{is2\pi j/N}.$$

$$⟨{\psi }_{\beta }^{\text{textsl}CFT}|{\tilde{\mathcal{O}}}^{\text{textsl}CFT,s}|0^{\text{textsl}CFT}⟩=\sum _{\alpha \in \mathcal{A}}{a}_{\alpha }⟨{\psi }_{\beta }^{\text{textsl}CFT}|{\psi }_{\alpha }^{\text{textsl}CFT,s}|0^{\text{textsl}CFT}⟩.$$

$$f_N^{\mathcal{O}}(\{a_{\alpha }\})$$

$$=$$

$$\mbox{(identity)}\leavevmode\nobreak\ \mathbb{1},\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \mbox{(stress tensor)}\leavevmode\nobreak\ T,\leavevmode\nobreak\ \bar{T}\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\$$

$$\displaystyle\mbox{(spin)}\leavevmode\nobreak\ \sigma,\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \partial_{\tau}\sigma,\leavevmode\nobreak\ \partial_{x}\sigma,\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \partial_{\tau}^{2}\sigma,\leavevmode\nobreak\ \partial_{x}\partial_{\tau}\sigma,\leavevmode\nobreak\ \partial_{x}^{2}\sigma,\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\$$

$$\displaystyle\mbox{(energy density)}\leavevmode\nobreak\ \epsilon,\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \partial_{\tau}\epsilon,\leavevmode\nobreak\ \partial_{x}\epsilon,\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \partial_{\tau}^{2}\epsilon,\leavevmode\nobreak\ \partial_{x}\partial_{\tau}\epsilon,\leavevmode\nobreak\ \partial_{x}^{2}\epsilon.\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\$$

$$\sigma^{\textsl{\tiny CFT}}\leavevmode\nobreak\ \stackrel{{\scriptstyle\tiny\mbox{approx}}}{{\sim}}\leavevmode\nobreak\ \mu X,\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\ \mu\approx 1.24514,\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\$$

$$\sigma^{\textsl{\tiny CFT}}\leavevmode\nobreak\ \stackrel{{\scriptstyle\tiny\mbox{approx}}}{{\sim}}\leavevmode\nobreak\ \mu^{\prime}X+\nu(XZ+ZX),\leavevmode\nobreak\ \leavevmode\nobreak\ \leavevmode\nobreak\$$

$$\text{leavevmode}\text{leavevmode}{ϵ}^{\text{textsl}CFT}\text{leavevmode}\stackrel{\text{mbox}approx}{\sim }\text{leavevmode}\mu (XX-Z),\text{leavevmode}\text{leavevmode}\text{leavevmode}\text{leavevmode}\mu \approx 1.5708,$$

$$\text{leavevmode}\text{leavevmode}{ϵ}^{\text{textsl}CFT}\text{leavevmode}\stackrel{\text{mbox}approx}{\sim }\text{leavevmode}{\mu }^{\prime }(XX-Z)+\nu (YY+XZX),$$

$$T^{\text{textsl}CFT},{\bar{T}}^{\text{textsl}CFT}\stackrel{\text{mbox}approx}{\sim }\mu (XX+Z)\pm \nu (XY-YX)+{\nu }^{\prime }1,$$

$$C_{\sigma \sigma ϵ}\approx 0.500000,\text{leavevmode}\text{leavevmode}{C}_{\sigma ϵ\sigma }\approx 0.50000,$$

$$X(j)\to -X(j),\text{leavevmode}\text{leavevmode}Y(j)\to -Y(j),\text{leavevmode}\text{leavevmode}Z(j)\to Z(j),$$

$$S_{j}X(j)\to S_{j}Y(j),\text{leavevmode}\text{leavevmode}\text{leavevmode}\text{leavevmode}{S}_{j}Y(j)\to S_{j+1}X(j+1),$$

$$z=e^{2\pi w/L},$$

$$\varphi (\tau ,x)$$

$$=$$

$$\varphi (z,\bar{z})=\sum _{m\in Z-h,\bar{m}\in Z-\bar{h}}{z}^{-m-h}{\bar{z}}^{-\bar{m}-\bar{h}}{\varphi }_{m,\bar{m}},$$

$$\varphi (0,x)={\left(\frac{2\pi }{L}\right)}^{\Delta }\sum _{m,\bar{m}\in Z}{e}^{2\pi i(m-\bar{m}+s_{\varphi })x/L}{\varphi }_{-m-h,-\bar{m}-\bar{h}},$$

Full text

Critical phenomena Cardy (1996); Sachdev (1999), characterized by universal behaviour and scale invariance, are of broad interest in various fields of physics, including statistical mechanics, condensed matter, and quantum fields. At a second order phase transition, scale invariance is often enhanced to a larger symmetry group, the conformal group, and the universal, low-energy physics of the critical system is then described by a conformal field theory (CFT) Polyakov (1970); Belavin et al. (1984); Friedan et al. (1984); Di Francesco et al. (2012), a field theory that is rigidly constrained by conformal invariance. This is particularly manifest in two dimensions, where the entire CFT is specified by a limited set of conformal data: its central charge $c^{\textsl{\tiny CFT}}$ together with the list of conformal dimensions $(h^{\textsl{\tiny CFT}}_{\alpha},\bar{h}^{\textsl{\tiny CFT}}_{\alpha})$ and operator product expansion (OPE) coefficients $C^{\textsl{\tiny CFT}}_{\alpha\beta\gamma}$ for its primary operators $\phi^{\textsl{\tiny CFT}}_{\alpha}$. From the conformal data one can derive the critical exponents of the phase transition as well as the rest of its universal properties.

There has been enormous progress in identifying possible conformal data for 1+1D CFTs, leading for instance to the famous characterization of unitary minimal models Belavin et al. (1984); Friedan et al. (1984). However, given a microscopic Hamiltonian $H=\sum_{j=1}^{N}h(j)$ for a critical quantum spin chain, numerically computing the conformal data of the emergent CFT remains a challenging task. At the core of the problem is the exponential growth of the Hilbert space with the size of the spin chain. Broadly speaking, two possible computational strategies are available. The first one is based on evaluating ground-state two-point and three-point correlators, which directly yield the conformal dimensions and the OPE coefficients, respectively Degli et al. (2004); Tagliacozzo et al. (2008); Xavier (2010); Evenbly and Vidal (2013); Stojevic et al. (2015); Evenbly and Vidal (2015). A second strategy, outlined by Cardy in the 80s Cardy (1984, 1986a), is based instead on exploiting the CFT operator-state correspondence Di Francesco et al. (2012). This correspondence, denoted $\psi_{\alpha}^{\textsl{\tiny CFT}}\sim|\psi_{\alpha}^{\textsl{\tiny CFT}}\rangle$, relates each scaling operator $\psi_{\alpha}^{\textsl{\tiny CFT}}$ of the CFT with a simultaneous eigenvector $|\psi_{\alpha}^{\textsl{\tiny CFT}}\rangle$ of the Hamiltonian and momentum operators $H^{\textsl{\tiny CFT}}$ and $P^{\textsl{\tiny CFT}}$ of the CFT on a circle of length $L$. In particular,

[TABLE]

where the energy and momentum of the state $|\psi^{\textsl{\tiny CFT}}_{\alpha}\rangle$ are given in terms of the scaling dimension $\Delta^{\textsl{\tiny CFT}}_{\alpha}=h^{\textsl{\tiny CFT}}_{\alpha}+\bar{h}^{\textsl{\tiny CFT}}_{\alpha}$ and conformal spin $s^{\textsl{\tiny CFT}}_{\alpha}=h^{\textsl{\tiny CFT}}_{\alpha}-\bar{h}^{\textsl{\tiny CFT}}_{\alpha}$ of the scaling operator $\psi_{\alpha}^{\textsl{\tiny CFT}}$. We can then exploit two classic results: (i) As first pointed out by Cardy Cardy (1984); Blöte et al. (1986); Cardy (1986b, a); Affleck (1986), the low energy states $|\psi_{\alpha}\rangle$ of a critical quantum spin chain on the circle are in one-to-one correspondence with CFT states, $|\psi_{\alpha}\rangle\sim|\psi^{\textsl{\tiny CFT}}_{\alpha}\rangle$, and approximately reproduce the spectrum of energies and momenta (1)-(2), from which we can estimate $c^{\textsl{\tiny CFT}}$ and $\Delta_{\alpha}^{\textsl{\tiny CFT}}$ and $s_{\alpha}^{\textsl{\tiny CFT}}$ (or $h_{\alpha}^{\textsl{\tiny CFT}}$ and $\bar{h}_{\alpha}^{\textsl{\tiny CFT}}$); (ii) Following Koo and Saleur Koo and Saleur (1994); Read and Saleur (2007); Dubail et al. (2010); Vasseur et al. (2012); Gainutdinov et al. (2013); Gainutdinov and Vasseur (2013); Bondesan et al. (2015), a Fourier expansion of the lattice Hamiltonian term $h(j)$ results in an approximate lattice realization of the CFT Virasoro generators. These can then be used Milsted and Vidal (2017) to identify the specific low energy states of the spin chain, denoted $|\phi_{\alpha}\rangle$, that correspond to CFT primary states $|\phi^{\textsl{\tiny CFT}}_{\alpha}\rangle$, which are the ones contributing to the conformal data. In addition, in order to significantly reduce finite size errors in the resulting numerical estimates of $c^{\textsl{\tiny CFT}}$ and $h^{\textsl{\tiny CFT}}_{\alpha},\bar{h}^{\textsl{\tiny CFT}}_{\alpha}$, in Ref. Zou et al. (2018) we recently demonstrated the use of periodic uniform matrix product states (puMPS) Rommer and Östlund (1997); Pirvu et al. (2012) to study critical chains made of up to several hundreds of quantum spins. Overall this second strategy, based on the operator-state correspondence, produces more systematic and accurate results than the approaches Degli et al. (2004); Tagliacozzo et al. (2008); Xavier (2010); Evenbly and Vidal (2013); Stojevic et al. (2015); Evenbly and Vidal (2015) based on two-point and three-point correlators. However, its major drawback is that its does not yield the OPE coefficients $C_{\alpha\beta\gamma}^{\textsl{\tiny CFT}}$.

In this paper we explain how to identify each local lattice operator $\mathcal{O}$, acting on the spin chain, with a corresponding linear combination of CFT scaling operators,

[TABLE]

More specifically, we show how to numerically compute the first few dominant terms in this expansion, corresponding to the CFT operators with the smallest scaling dimensions. As a main application, we then explain how to extract a lattice estimate $C_{\alpha\beta\gamma}$ of the OPE coefficients $C^{\textsl{\tiny CFT}}_{\alpha\beta\gamma}$, by computing the matrix element

[TABLE]

where $\phi_{\alpha}$ is an approximate lattice realization of the CFT primary operator $\phi^{\textsl{\tiny CFT}}_{\alpha}$ and $|\phi_{\beta}\rangle$ and $|\phi_{\gamma}\rangle$ are a pair of primary states of the critical spin chain. In this way we successfully complete Cardy’s ambitious program to extract conformal data from a critical lattice Hamiltonian $H$ by exploiting the operator-state correspondence. We demonstrate the approach, valid for any quantum spin chain, by computing the leading terms of the expansion (3) for all one-site and two-site operators of the critical Ising model, see Table 1, as well as its non-trivial OPE coefficient $C^{\textsl{\tiny CFT}}_{\sigma\sigma\epsilon}$. We also briefly enumerate other future applications of the method.

Exciting the CFT vacuum with local operators.— We start by reviewing some basic facts. A $1+1$ dimensional CFT can be quantized on a cylinder $S^{1}\times\mathbb{R}$, where the compactified dimension represents space, with coordinate $x\in[0,L)$, and the other dimension represents Euclidean time, with coordinate $\tau\in\mathbb{R}$. On the $\tau=0$ circle we build the Hilbert space, spanned by the states $|\psi_{\alpha}^{\textsl{\tiny CFT}}\rangle$ in Eqs. (1)-(2). Let $\mathcal{O}^{\textsl{\tiny CFT}}(x)$ denote a local operator acting on this circle, with Fourier mode decomposition

[TABLE]

Applying Fourier mode $\mathcal{O}^{\textsl{\tiny CFT},s}$ on the vacuum $|0^{\textsl{\tiny CFT}}\rangle$ results in an eigenstate of $P^{\textsl{\tiny CFT}}$ with momentum $(2\pi/L)s$. When $\mathcal{O}^{\textsl{\tiny CFT}}(x)$ is a primary field $\phi^{\textsl{\tiny CFT}}_{\alpha}(x)$, we obtain Sup

[TABLE]

where $\phi^{\textsl{\tiny CFT}}_{\alpha,(m,\bar{m})}(x)\equiv\partial^{m}\bar{\partial}^{\bar{m}}\phi^{\textsl{\tiny CFT}}_{\alpha}(x)$ denotes a derivative descendant with spin $s_{\alpha}^{\textsl{\tiny CFT}}+m-\bar{m}$, $\partial\equiv(\partial_{\tau}-i\partial_{x})/2$ is the complex derivative (with $\bar{\partial}=(\partial_{\tau}+i\partial_{x})/2$), and

[TABLE]

More generally, analogous expressions can be obtained when $\mathcal{O}^{\textsl{\tiny CFT}}$ is not a primary operator. Here we will use the specific case where $\mathcal{O}^{\textsl{\tiny CFT}}$ is a derivative descendant $\phi^{\textsl{\tiny CFT}}_{\alpha,(k,\bar{k})}$ of the primary $\phi^{\textsl{\tiny CFT}}_{\alpha}$, for which one finds Sup

[TABLE]

Finally, the OPE coefficients can be extracted from Sup

[TABLE]

where $\phi_{\alpha}^{\textsl{\tiny CFT}}$ is a primary operator and $|\phi_{\beta}^{\textsl{\tiny CFT}}\rangle$ and $|\phi_{\gamma}^{\textsl{\tiny CFT}}\rangle$ are a pair of CFT primary states on the circle.

Lattice operators as CFT scaling operators.— Consider now a critical quantum spin chain and a local operator $\mathcal{O}$, acting on a small number of spins, to which we would like to assign a linear combination of CFT scaling operators as in Eq. (3). In practice we will produce an approximate, truncated expansion of the form

[TABLE]

using only operators $\psi^{\textsl{\tiny CFT}}_{\alpha}$ in a preselected finite set $\mathcal{A}$. By optimizing the coefficients $a_{\alpha}$ (see below), we hope to obtain a truncated expansion (11) such that

[TABLE]

where the matrix elements are between low energy states $|\psi_{\alpha}\rangle$ and $|\psi_{\beta}\rangle$, we have equated the size $N$ of the spin chain with the size $L$ of the CFT circle, and $\Delta_{c}$ is the lowest scaling dimension among operators not included in $\mathcal{A}$. Thus, the accuracy of the expansion should systematically improve (the subleading finite-size corrections be further reduced) by adding more scaling operators in $\mathcal{A}$.

Constraining $\tilde{\mathcal{O}}^{\textsl{\tiny CFT}}$ without knowing the OPE.— In analogy with Eqs. (5)-(6), we first Fourier expand $\mathcal{O}$,

[TABLE]

Given a finite set $\mathcal{B}$ of low energy states $|\psi_{\beta}\rangle$ and a range $\mathcal{S}$ of values $s$, we can numerically evaluate the matrix elements $b_{\beta,s}\equiv\langle\psi_{\beta}|\mathcal{O}^{s}|0\rangle$ between the spin chain ground state $|0\rangle$ and state $|\psi_{\beta}\rangle$, for all $|\psi_{\beta}\rangle\in\mathcal{B}$ and $s\in\mathcal{S}$. With the ability to analytically compute $B^{\alpha}_{\beta,s}\equiv\langle\psi_{\beta}^{\textsl{\tiny CFT}}|\psi_{\alpha}^{\textsl{\tiny CFT},s}|0^{\textsl{\tiny CFT}}\rangle$ (using e.g. Eqs. (7)-(9)), we can also evaluate the corresponding CFT matrix elements

[TABLE]

In this way we can search for the coefficients $a_{\alpha}$ such that $\langle\psi_{\beta}^{\textsl{\tiny CFT}}|\tilde{\mathcal{O}}^{\textsl{\tiny CFT},s}|0^{\textsl{\tiny CFT}}\rangle$ best approximates $\langle\psi_{\beta}|\tilde{\mathcal{O}}^{s}|0\rangle$ for all relevant $s$ and $\beta$, by minimizing the cost function

[TABLE]

Importantly, $f_{N}^{\mathcal{O}}(\{a_{\alpha}\})$ depends only on matrix elements involving the vacuum and one excited state (analogous to a CFT two-point correlator) and not on matrix elements involving two excited states (analogous to a three-point correlator), so that it does not require any knowledge of the OPE coefficients $C_{\alpha\beta\gamma}^{\textsl{\tiny CFT}}$.

Optimization.— The algorithm is divided into parts I and II, involving states and operators, respectively.

Part I (low energy states): taking the critical Hamiltonian $H=\sum_{j=1}^{N}h(j)$ of a periodic spin chain as the only input, we compute low energy eigenstates of $H$, using e.g. exact diagonalization for small $N$ and puMPS for larger $N$ Zou et al. (2018). We then use the techniques of Refs. Milsted and Vidal (2017); Zou et al. (2018) to (i) for each low energy state $|\psi_{\alpha}\rangle\sim|\psi_{\alpha}^{\textsl{\tiny CFT}}\rangle$, obtain estimates $\Delta_{\alpha}$ and $s_{\alpha}$ for its scaling dimension and conformal spin; (ii) identify primary states $|\phi_{\alpha}\rangle\sim|\phi^{\textsl{\tiny CFT}}_{\alpha}\rangle$ and their descendants, thus organizing the low energy states into conformal towers. The above tasks involve a large-$N$ extrapolation. At this point we make a judicious choice of sets $\mathcal{A}$, $\mathcal{B}$, and $\mathcal{S}$ (see example below).

Part II (operators): For a fixed system size $N$, we compute the CFT matrix elements $\langle\psi_{\beta}^{\textsl{\tiny CFT}}|\psi_{\alpha}^{\textsl{\tiny CFT},s}|0^{\textsl{\tiny CFT}}\rangle$ for all $\psi_{\alpha}^{\textsl{\tiny CFT}}\in\mathcal{A}$, $|\psi_{\beta}^{\textsl{\tiny CFT}}\rangle\in\mathcal{B}$ and $s\in\mathcal{S}$ using the corresponding analytical expressions. However, here we employ the previously estimated conformal dimensions $h_{\alpha}\equiv(\Delta_{\alpha}+s_{\alpha})/2,\bar{h}_{\alpha}\equiv(\Delta_{\alpha}-s_{\alpha})/2$ instead of their unknown exact values $h_{\alpha}^{\textsl{\tiny CFT}},\bar{h}_{\alpha}^{\textsl{\tiny CFT}}$, so that no previous knowledge of the emergent CFT is required. Then, for each choice of lattice operator $\mathcal{O}$, we compute $\langle\psi_{\beta}|\mathcal{O}^{s}|0\rangle$ for all $\psi_{\beta}\in\mathcal{B}$ and $s\in\mathcal{S}$, and use linear least-squares regression to minimize $f_{N}^{\mathcal{O}}(\{a_{\alpha}\})$ in Eq. (15), resulting in a set of optimal coefficients $a_{\alpha}(N)$. Finally, we repeat the entire calculation for several values of $N$ and extrapolate to $N\rightarrow\infty$, which results in the coefficients $a_{\alpha}$ in (11).

Example: Critical quantum Ising model.— To illustrate the approach we consider the critical transverse field Ising Hamiltonian, with local Hamiltonian term $h(j)=-X(j)X(j+1)-Z(j)$. We used puMPS with bond dimension in the range $12\leq D\leq 22$ to address systems of size $18\leq N\leq 48$, with the largest system size requiring several minutes on a laptop with 4 CPU (2.8 GHz) and 2 GB RAM. Part I above yields three conformal towers Zou et al. (2018), including the identity (or ground state) tower, present in all CFTs, and two additional conformal towers. We call the corresponding primaries spin $\sigma$ and energy density $\epsilon$ to follow the standard Ising CFT nomenclature (we reiterate, however, that no previous knowledge of the emergent CFT was used). For the set $\mathcal{A}$ we choose operators

[TABLE]

in the identity conformal tower (notice that $\mathbb{1}$, $T$ and $\bar{T}$ are present in any CFT) and the primaries and first and second derivative descendants in the other two towers,

[TABLE]

For $\mathcal{B}$ we choose all states $|\psi_{\beta}\rangle$ with scaling dimension $\Delta_{\alpha}\leq 3+1/8$, namely the 23 lowest energy states. After normalizing $H$ such that $\Delta_{\mathbb{1}}=0$ and $\Delta_{T}=\Delta_{\bar{T}}=2$ Milsted and Vidal (2017), our estimates for the scaling dimensions of primaries are $\Delta_{\sigma}=0.1249995$ and $\Delta_{\epsilon}=0.9999994$, together with exact values $s_{\mathbb{1}}=s_{\sigma}=s_{\epsilon}=0$, $s_{T}=-s_{\bar{T}}=2$ for the conformal spins. Scaling dimensions and conformal spins of derivative descendants are obtained by adding integers to these values. Finally, for $\mathcal{S}$ we choose $-3\leq s\leq 3$.

In part II, we evaluated all matrix elements $\langle\psi^{\textsl{\tiny CFT}}_{\beta}|\psi_{\alpha}^{\textsl{\tiny CFT}}|0^{\textsl{\tiny CFT}}\rangle$ using both Eqs. (7)-(9) and the scaling dimensions and conformal spins quoted above. Then we optimized the truncated expansion (11) for each Pauli operator $X,Y,Z$ acting on a single spin, for pairs of Pauli operators acting on two continuous spins, etc, see Table 1 and Fig. 1. Some of the coefficients $a_{\alpha}$ reproduce up to 6 significant digits of their exact value, obtained using the free fermion representation of the Ising model Sup .

Primary operators and OPE coefficients.— By inverting the relations in Table 1, we can build linear combinations of lattice operators whose leading contribution is a targeted CFT scaling operator $\psi^{\textsl{\tiny CFT}}_{\alpha}$. For instance, if the target is the spin primary $\sigma^{\textsl{\tiny CFT}}$, its simplest realization is given in terms of the Pauli matrix $X$

[TABLE]

where the approximation can be seen from Table 1 to include both $\partial_{x}^{2}\sigma^{\textsl{\tiny CFT}}$ and $\partial_{\tau}^{2}\sigma^{\textsl{\tiny CFT}}$ as subleading corrections. An improved lattice realization is then given by

[TABLE]

where $\nu\approx-0.026$, $\mu^{\prime}\approx 1.27$ and, importantly, the subleading corrections due to $\partial_{\tau}^{2}\sigma^{\textsl{\tiny CFT}}$ have been eliminated This is particularly useful below for the computation of $C_{\sigma\sigma\epsilon}$, seen to be insensitive to the correction $\partial_{x}^{2}\sigma^{\textsl{\tiny CFT}}$ still present in (21). Similarly, for the primary $\epsilon^{\textsl{\tiny CFT}}$ we find

[TABLE]

with $\partial_{\tau}^{2}\epsilon^{\textsl{\tiny CFT}}$ as a subleading correction, and the improved

[TABLE]

for $\nu\approx-0.19$ and $\mu^{\prime}\approx 1.76$, with no subleading corrections in $\mathcal{A}$. Finally, we can also obtain the stress tensor

[TABLE]

where $\mu\approx-1.571$, $\nu\approx 0.7854$, $\nu^{\prime}\approx 2.000$ and, again, there are no subleading contributions in $\mathcal{A}$.

Equipped with a lattice realization of the primary operators $\sigma^{\textsl{\tiny CFT}}$ and $\epsilon^{\textsl{\tiny CFT}}$, we can finally use Eq. (4) (motivated by Eq. (10)) to compute estimates $C_{\alpha\beta\gamma}$ for the OPE coefficients $C^{\textsl{\tiny CFT}}_{\alpha\beta\gamma}$. The non-zero coefficients are

[TABLE]

obtained by computing $\langle\sigma|\mathcal{O}_{\sigma}(0)|\epsilon\rangle$ and $\langle\sigma|\mathcal{O}_{\epsilon}(0)|\sigma\rangle$, with $\mathcal{O}_{\sigma}$ and $\mathcal{O}_{\epsilon}$ the lattice operators in Eqs. (21) and (23) and extrapolating to the thermodynamic limit Sup . Fig. 2 shows a change in scaling from $N^{-2}$ to $N^{-4}$ in the error of $C_{\sigma\sigma\epsilon}$ when replacing the lattice realization (20) of $\sigma^{\textsl{\tiny CFT}}$ with the improved lattice realization (21).

Symmetry, duality, and free fermions.— The critical Ising model is invariant under a global spin flip

[TABLE]

and (on the real line) under the Kramers-Wannier duality

[TABLE]

where $S_{j}\equiv\prod_{i<j}Z(i)$, such that e.g. $X(j)X(j+1)\rightarrow Z(j+1)$ and $Z(j)\rightarrow X(j)X(j+1)$. Our computations above, which made no use of these properties, yield results fully consistent with them, as should be expected. For instance, all CFT and lattice operators in Eqs. (20)-(21) are odd under spin flip symmetry; operators in Eqs. (22)-(23) are even under spin flip and odd under the duality transformation; operators in Eq. (24) are even under spin flip and duality transformation; finally, the OPE coefficients (25) are the only non-vanishing ones allowed by symmetry. Moreover, the spin model admits a free fermion representation. This allows for an exact computation of some of the coefficients in Table 1 Sup , which we used to confirm the accuracy of the numerical results.

Discussion.— Given a critical quantum spin chain Hamiltonian $H$ as the only input, in this paper we have explained how to identify a local lattice operator with a corresponding expansion (11) in terms of scaling operators $\psi^{\textsl{\tiny CFT}}_{\alpha}$ of the emergent CFT. As demonstrated for the critical Ising model, this allows us to build lattice versions of specific CFT scaling operators. In particular, by targeting primary operators, one can compute OPE coefficients, thereby completing Cardy’s program to numerically extract the conformal data from low-energy states of a critical spin chain by exploiting the operator-state correspondence. Our approach, which can be extended to address non-local scaling operators Zou et al. (2019), has other useful applications, explored in subsequent work (see also Sup ). For instance, in a generic critical quantum spin chain one can now modify the original Hamiltonian by adding relevant (or irrelevant) scaling operators on demand. Then, using e.g. the techniques demonstrated in Ref. Zou et al. (2018), one can study, fully non-perturbatively, the renormalization group flow away from (respectively, back to) the initial CFT, along pre-determined directions. Conversely, given a near-critical lattice Hamiltonian, one can tune it closer to criticality by removing relevant perturbations from it.

Acknowlegments — We thank Qi Hu and Martin Ganahl for many useful discussions, as well as Alexander Zamolodchikov for valuable comments. The authors acknowledge support from the Simons Foundation (Many Electron Collaboration) and Compute Canada. Research at Perimeter Institute is supported by the Government of Canada through the Department of Innovation, Science and Economic Development Canada and by the Province of Ontario through the Ministry of Research, Innovation and Science.

Appendix A General CFT results

In this section, we derive CFT matrix elements that are used in the main text. While in the main text we have used the superscript ${}^{\textsl{\tiny CFT}}$ to indicate that a state or operator belong to the CFT, in this section we omit the superscript if it is clear from the context that we are referring to a CFT quantity.

A.1 CFT matrix elements involving the ground state

We are first concerned with matrix elements involving the ground state and one excited state. They are equivalent to two-point correlations and do not depend on the OPE coefficients of the CFT.

A 1+1 dimensional CFT can be quantized on a cylinder with complex coordinates $w=\tau+ix$ and $\bar{w}=\tau-ix$, where $-\infty<\tau<+\infty$ represents the imaginary time coordinate and $0\leq x<L$ represents the spatial (periodic) coordinate. A primary field $\phi(\tau,x)$ can be mapped to the complex plane with coordinates $(z,\bar{z})$ by the conformal transformation

[TABLE]

which using that $dz/dw=2\pi z/L$ leads to

[TABLE]

where $h,\bar{h}$ are the conformal dimensions of $\phi(z,\bar{z})$, and $\Delta=h+\bar{h}$ is its scaling dimension.

On the complex plane, the field can be Laurent expanded around the origin,

[TABLE]

which becomes Fourier expansion on a time slice of the cylinder,

[TABLE]

where $s_{\phi}=h-\bar{h}$ is the conformal spin.

Applying a Fourier mode on the ground state gives

[TABLE]

Therefore, acting with $\phi(0,x)$ on the ground state creates the state $|\phi\rangle$ and all its derivative descendants,

[TABLE]

where $|\phi_{(m,\bar{m})}\rangle$ is a normalized state with conformal dimensions $(m+h,\bar{m}+\bar{h})$ and

[TABLE]

is a normalization constant, which can be obtained by successively applying the Virasoro algebra.

We can then obtain the matrix elements of the Fourier mode

[TABLE]

by combining Eqs. (32)-(36). Adding the superscript ${}^{\textsl{\tiny CFT}}$ and the subscript α to the fields we finally obtain the equation

[TABLE]

where

[TABLE]

as used in the main text.

Let us generalize the above equations to Fourier modes of descendant operators. Note that the transformation rule Eq. (29) does not apply. Descendant fields are defined through OPE

[TABLE]

where in particular, for $n=-2$ and $n=-1$ we have singular terms that involve $L_{0}\phi=h\phi$ and $L_{-1}\phi=\partial\phi$. For concreteness, we next consider the simplest descendant $L_{-1}\phi$. We can extract it from the OPE by a contour integral

[TABLE]

Note that the operator product should be understood as a time-ordered product. We can deform the coutour to two closed circles $\tau_{\pm}=\mathrm{Re}w^{\prime}\pm\epsilon$ to obtain $\partial\phi(\tau^{\prime},x^{\prime})=$

[TABLE]

By transforming $T(z)$ onto the cylinder (analogous to Eq. (29) with another Schwarzian term proportional to the central charge on the RHS), we obtain

[TABLE]

for any $\tau$. Thus we arrive at

[TABLE]

Taking the Fourier mode $s$ on both sides, setting $\tau=0$, and acting on the ground state, we obtain an expression similar to Eq. (37) (we also add superscripts and subscripts to compare with the main text):

[TABLE]

and more generally

[TABLE]

For later use, let us consider a special case $\partial_{x}\phi_{\alpha}^{\textsl{\tiny CFT}}\equiv i(\partial-\bar{\partial})\phi_{\alpha}^{\textsl{\tiny CFT}}$,

[TABLE]

where the first equality follows from linear combing the cases of $k=1,\bar{k}=0$ and $k=0,\bar{k}=1$ in Eq. (45) and $s^{\textsl{\tiny CFT}}_{\alpha}=h^{\textsl{\tiny CFT}}_{\alpha}-\bar{h}^{\textsl{\tiny CFT}}_{\alpha}$, and the second equality follows from $\delta_{m-\bar{m}+s^{\textsl{\tiny CFT}}_{\alpha},s}$ in Eq. (38).

Generic descendant operators are composite operators of the stress tensors and the primary operators. They can create states (quasi-primary states) which are not derivative descendant states. An analogous formula can still be derived, but is far more complicated. However, the case for the stress tensor scaling operators $T^{\textsl{\tiny CFT}}$ and $\bar{T}^{\textsl{\tiny CFT}}$ (quasi-primary operators in the conformal tower of the identity) can be derived in a simpler way, since their Fourier modes are Virasoro generators,

[TABLE]

Applying the $L^{\textsl{\tiny CFT}}_{-s},\bar{L}^{\textsl{\tiny CFT}}_{-s}$ to the ground state produces

[TABLE]

for $s\geq 2$, where

[TABLE]

is the normalization constant. Combining the equations above gives the matrix elements of $T^{\textsl{\tiny CFT},s},\bar{T}^{\textsl{\tiny CFT},s}$.

In the case of the Ising CFT (see next appendix), all descendants of $\sigma^{\textsl{\tiny CFT}}$ and $\epsilon^{\textsl{\tiny CFT}}$ below level $2$ are derivative descendants. Together with the stress tensor operators, these are all the scaling operators included in the truncated set $\mathcal{A}$. Therefore, the above equations are all that is needed for the particular application considered in the main text.

A.2 OPE coefficients

On the complex plane, the OPE coefficients are defined by

[TABLE]

Transforming $\phi^{\textsl{\tiny CFT}}_{\beta}$ onto the cylinder with Eq. (29), we have

[TABLE]

Since CFT states on the circle $\tau=0$ of the cylinder are eigenstates of the translation $x\rightarrow x+\delta x$, only the Fourier mode of $\phi^{\textsl{\tiny CFT}}_{\beta}$ with momentum $s_{\alpha}-s_{\gamma}$ contributes to the above equation. Therefore,

[TABLE]

Note that generally the OPE coefficient is invariant under even permutations of its three labels (e.g. $C_{\alpha\beta\gamma}=C_{\beta\gamma\alpha}$), and becomes complex conjugated under odd permutations (e.g. $C_{\alpha\beta\gamma}=C_{\beta\alpha\gamma}^{*}$).

Any other three point function of the same operators is related to the standard form Eq. (53) by some conformal transformation. Therefore, $C^{\textsl{\tiny CFT}}_{\alpha\beta\gamma}$ determines all three point correlation functions. For example,

[TABLE]

which follows from Eq. (43) and its antiholomorphic analogue. Matrix elements of general descendants can also be derived in a more complicated way, which we omit here since we will not need them for the particular application in the main text.

Appendix B The Ising CFT

The Ising CFT describes the low energy, long distance, universal physics of many critical lattice models, including e.g. the critical Ising quantum spin chain and the quantum critical axial next-to-nearest neighbor Ising quantum spin chain. The Ising CFT is also the first one in the series of unitary minimal models, which can be solved exactly. In this section we shall review some properties of the Ising CFT that are used in this paper. Unless otherwise stated, objects in this section are CFT objects and we omit the superscript ${}^{\textsl{\tiny CFT}}$.

B.1 Conformal data

As any other unitary minimal models, this CFT has a finite number of primary operators, resulting in a finite amount of conformal data. Specifically, the central charge is $c=1/2$ and there are three primary fields, namely the identity operator $\mathbb{1}$ (present in any CFT), the spin operator $\sigma$ and the energy density operator $\sigma$. They all have conformal spin $s_{\mathbb{1}}=s_{\sigma}=s_{\epsilon}=0$ and their scaling dimensions are $\Delta_{\mathbb{1}}=0$, $\Delta_{\sigma}=1/8$, and $\Delta_{\epsilon}=1$, respectively. The only nonzero OPE coefficients (up to permutations of the indices) are $C_{\alpha\beta\mathbb{1}}=\delta_{\alpha\beta}$ and $C_{\sigma\sigma\epsilon}=1/2$.

B.2 Null fields and conformal towers

The Ising CFT has null fields, whose correlation functions are zero and therefore they do not correspond to a state in the CFT. These are

[TABLE]

As a result, all descendants of $\sigma$ and $\epsilon$ under level $2$ are derivative descendants. The lowest descendant that is not a derivative descendant is $L_{-3}\sigma$ in the $\sigma$ tower and $L_{-4}\epsilon$ in the $\epsilon$ tower (and those with $\bar{L}_{-n}$). We shall see later that the $L_{-4}\epsilon$ operator is responsible for the finite-size corrections of the OPE coefficient $C_{\sigma\epsilon\sigma}$ computed on the lattice.

Appendix C General critical lattice models

In this section we consider critical lattice models. Given an operator on the lattice, we investigate how to identify it with a sum of CFT operators in the continuum.

C.1 Fourier modes of multi-site operators

Recall that the Fourier mode $\mathcal{O}^{s}$ of a lattice operator $\mathcal{O}$ is defined by

[TABLE]

For a one-site operator at site $j$, the position assignment $x_{j}=j$ appears uncontroversial. However, for an operator $\mathcal{O}(j)$ supported on multiple sites, say from site $j$ to site $j+n$, the position $x_{j}$ is not uniquely determined. We have to decide how to assign a specific position $x_{j}\in(j,j+n)$ to it. Different assignments will lead to different expansions in terms of CFT operators. However, it can be shown that any two such expansions have the same dominant CFT scaling operator, and the difference between the two expansions is dominated by the derivative of this dominant CFT operator.

Let us illustrate the above with an example for the critical Ising model. Consider $\mathcal{O}_{1}(j)=-X(j)X(j+1)$. We have seen that its CFT expansion $\mathcal{O}^{\textsl{\tiny CFT}}_{1}$ includes both $\mathbb{1}^{\textsl{\tiny CFT}}$ and $\epsilon^{\textsl{\tiny CFT}}$ contributions. For this lattice operator, we may assign e.g. $x_{j}=j$, $x_{j}=j+1/2$, pr $x_{j}=j+1$. Only the second choice preserves spatial parity, and therefore no $\partial_{x}\epsilon^{\textsl{\tiny CFT}}$ term is allowed in the expansion of $O^{\textsl{\tiny CFT}}_{1}$. The other two choices would result in a $\partial_{x}\epsilon^{\textsl{\tiny CFT}}$ term in the expansion of $\mathcal{O}^{\textsl{\tiny CFT}}_{1}$, which is in accordance with the fact that our assignment of position has explicitly broken spatial parity. Nevertheless, the expansion coefficients in front of $\mathbb{1}^{\textsl{\tiny CFT}}$ and $\epsilon^{\textsl{\tiny CFT}}$ are independent of our assignment of position $x_{j}$.

The specific assignment $x_{j}=j+1/2$ for $\mathcal{O}_{1}(j)$ is important when combining $\mathcal{O}_{1}(j)=-X(j)X(j+1)$ with $\mathcal{O}_{2}(j)=-Z(j)$ to form the Hamiltonian density $h(j)=\mathcal{O}_{1}(j)+\mathcal{O}_{2}(j)$. In order for the Fourier mode $h^{s}$ to correspond to a linear combination of Virasoro generators $L^{\textsl{\tiny CFT}}_{-s}+\bar{L}^{\textsl{\tiny CFT}}_{s}$, it has been shown numerically Milsted and Vidal (2017) that the correct choice is $x_{j}=j+1/2$ for $O_{1}(j)$ and $x_{j}=j$ for $O_{2}(j)$. If we have chosen a different $x_{j}$ for $O_{1}(j)$, the Fourier mode $h^{s}$ $(s\neq 0)$ would connect states in identity tower and $\epsilon$ tower. This is exactly the consequence of the $\partial_{x}\epsilon^{\textsl{\tiny CFT}}$ term in the expansion of $h$.

In the main text, for $Z_{2}$ odd operators we have used a simple middle point rule to assign positions, that is, an operator with support from site $j$ to $j+n$ is assigned position $x_{j}=j+n/2$. $Z_{2}$ even operators are assigned positions at the middle point of their fermionic representations in Table 2. The fermions $\psi(2j-1)$ and $\psi(2j)$ are assigned positions $x=j-1/4$ and $x=j+1/4$ respectively (according to the picture that one spin degree of freedom splits into two Majorana fermion degree of freedom). For most $Z_{2}$ even operators that are considered here, it coincides with the middle point rule in the spin representation. The exceptions are: $X(j)Y(j+1)$ is assigned position $x_{j}=j+3/4$ and $Y(j)X(j+1)$ is assigned position $x_{j}=j+1/4$.

We point out that different position assignments can also be chosen which are equally valid, as long as a consistent choice of convention is kept throughout the computations.

C.2 Fixing phases of low energy eigenstates

Diagonalization of the Hamiltonian yields a set of eigenstates with arbitrary complex phases. However, in the CFT calculations needed in our cost function (used in the main text to identify lattice operators with CFT operators), we are comparing lattice matrix elements with CFT matrix elements directly. To do this in a meaningful way, we first have to fix the complex phases of the low energy eigenstates of the critical spin chain using the same conventions used in the CFT. This is achieved by requiring certain matrix elements of lattice operators to have the same phases as in their CFT counterparts.

First, Fourier modes $h^{s}\sim L^{\textsl{\tiny CFT}}_{-s}+\bar{L}^{\textsl{\tiny CFT}}_{s}$ $(s\neq 0)$ of the lattice Hamiltonian density $h(j)$ are ladder operators in the scaling limit. In a CFT,

[TABLE]

Accordingly, we will require that the equivalent lattice matrix elements also satisfy

[TABLE]

up to finite-size corrections. In practice, for a given conformal tower, we fix relative phases between descendant states and the primary state $|\phi\rangle$ level by level. Starting with $|\psi_{\alpha}\rangle=|\phi\rangle$, we require the above matrix elements with $n=1$ and $n=2$ to be real and positive. Then we continue to $|\psi_{\alpha}\rangle=|L_{-n}\phi\rangle$ and $|\bar{L}_{-n}\phi\rangle$ ($n=1,2$) and fix the phases of higher level descendants. This is done until all the selected states in the cost function have their phases fixed with respect to the primary states.

In the remaining we would like to fix relative phases between primary states $|\phi_{\alpha}\rangle$. In the CFT,

[TABLE]

On the lattice, we first find an operator $O$ which has $\phi^{\textsl{\tiny CFT}}_{\alpha}$ in its expansion, and then require

[TABLE]

to be real and positive.

In the Ising CFT, primary fields are Hermitian,

[TABLE]

Therefore, for the critical Ising model, we can choose

[TABLE]

although we could have chosen other operators (for example, $XZ+ZX$ for $\sigma^{\textsl{\tiny CFT}}$ and $-Z$ for $\epsilon^{\textsl{\tiny CFT}}$, which turns out to be completely equivalent). After that, we have fixed the complex phases in all low energy eigenstates relative to the ground state.

C.3 OPE coefficients

Following the CFT expression, OPE coefficients can be computed on the lattice by

[TABLE]

where $\phi_{\beta}$ is the lattice operator corresponding to $\phi^{\textsl{\tiny CFT}}_{\beta}$. Since each state is an eigenstate of the translation operator, the above equation can also be written as

[TABLE]

because only the Fourier mode with momentum $s_{\alpha}-s_{\gamma}$ is consistent with momentum conservation.

Note that the exact lattice representation $\phi_{\beta}$ of a CFT scaling operator $\phi_{\beta}^{\textsl{\tiny CFT}}$ is not expected to be supported on a finite number of lattice sites in general. However, here we always work with a truncated set $\mathcal{A}$ of scaling operators and the resulting lattice operator $\mathcal{O}_{\phi_{\beta}}$ is only approximate but with finite local support and such that the above matrix elements differ by finite-size corrections.

After fixing the complex phases of lattice eigenstates as discussed above, for the Ising model we numerically find that

[TABLE]

and

[TABLE]

are both real and positive for all approximate lattice representations of primary operators that were used.

C.4 Sources of errors

In the main text, we have described the cost function that is minimized to obtain the truncated CFT operator expansion $\tilde{\mathcal{O}}^{\textsl{\tiny CFT}}$ corresponding to a lattice operator $\mathcal{O}$. Here we consider possible sources of errors in the construction, and how we can reduce the error in practice. There are three sources of errors.

The CFT operator space is truncated to a finite set $\mathcal{A}$. This precludes an exact correspondence between CFT operators and lattice operators. 2. 2.

The lattice has a finite number $N$ of sites, which leads to errors (due to subleading finite-size corrections) in the numerical estimates of the scaling dimensions and central charge used in order to evaluate CFT matrix elements in the cost function. 3. 3.

The numerical diagonalization of the Hamiltonian (e.g. using matrix product states MPS) produces approximate eigenstates (e.g. due to the finite bond dimension of the MPS).

It will be argued in this section that the first source causes errors in the expansion coefficients that decay as $1/N^{p}$, with the power depending on the truncated space $\mathcal{A}$ of CFT operators. The power-law convergence of expansion coefficients is later confirmed numerically by the results obtained with the Ising model. We also briefly comment on the other two sources of error, which are assumed to not be the dominant ones.

C.4.1 Errors due to a truncated set $\mathcal{A}$ of CFT operators

Given a truncated set $\mathcal{A}$ of CFT scaling operators $\{\psi^{\textsl{\tiny CFT}}_{\alpha}(x)\}$ and a lattice operator $\mathcal{O}$, we find the coefficients $a_{\alpha}(N)$ such that

[TABLE]

minimizes the cost function

[TABLE]

The exactly correspondent CFT operator $O^{\textsl{\tiny CFT}}$, which typically involves an infinite sum of scaling operators, satisfies

[TABLE]

for any $\beta$ and $s$. The goal is to estimate how far the coefficient $a_{\alpha}(N)$ that minimizes the cost function is away from $a_{\alpha}$.

Denote by $\mathcal{A}_{c}$ the set of scaling operators $\{\psi^{\textsl{\tiny CFT}}_{c,\alpha^{\prime}}(x)\}$ (c is a label for ”complementary”) that together with $\{\psi^{\textsl{\tiny CFT}}_{\alpha}(x)\}$ form a basis that expands $O^{\textsl{\tiny CFT}}(x)$, then we can expand

[TABLE]

Denote the difference $\delta a_{\alpha}(N)=a_{\alpha}-a_{\alpha}(N)$. Using Eqs. (73)-(75), we can express the cost function solely in terms of CFT matrix elements.

[TABLE]

For simplicity, we first consider the case where the expansion Eq. (76) only involves operators in one conformal tower. In the limit of large $N$, the second term scales as $N^{-\Delta_{c}}$, where $\Delta_{c}$ denotes the smallest scaling dimension of the operators in $\{\psi^{\textsl{\tiny CFT}}_{c,\alpha^{\prime}}(x)\}$ that have nonzero coefficient $b_{\alpha^{\prime}_{0}}\neq 0$. Moreover, the leading contribution of the second term cannot be completely eliminated by fine tuning $\delta a_{\alpha}(N)$, since otherwise $\psi^{\textsl{\tiny CFT}}_{c,\alpha^{\prime}_{0}}(x)$ would be a linear combination of $\{\psi^{\textsl{\tiny CFT}}_{\alpha}(x)\}$. Therefore, the minimum of the cost function $f^{O}_{N}$ scales as $N^{-2\Delta_{c}}$. Minimizing the cost function by fine-tuning $\delta a_{\alpha}(N)$ thus yields

[TABLE]

where $\Delta_{\alpha}$ is the scaling dimension of $\psi^{\textsl{\tiny CFT}}_{\alpha}$. In practice, we include in $\mathcal{A}$ all possible operators in $\{\psi^{\textsl{\tiny CFT}}_{\alpha}(x)\}$ up to scaling dimension $\Delta_{\mathrm{max}}$. Then by definition $\Delta_{c}>\Delta_{\mathrm{max}}$. Therefore, the error becomes smaller as we include operators in $\mathcal{A}$ with higher scaling dimensions so as to increase $\Delta_{\mathrm{max}}$. Another way to reduce error is to go to large sizes, if this error is dominant over other sources of error mentioned below.

If the expansion Eq. (76) involves operators in different conformal towers, then in Eq. (77) the sum over $\beta$ splits into different conformal towers. For each conformal tower, the sum over $\alpha$ and $\alpha^{\prime}$ are restricted to the same conformal tower. Following the same arguments, we define $\Delta_{c}$ for each conformal tower as the smallest scaling dimension in $\{\psi^{\textsl{\tiny CFT}}_{c,\alpha^{\prime}}(x)\}$ in that conformal tower, and Eq. (78) still holds for operators $\psi^{\textsl{\tiny CFT}}_{\alpha}$ in that conformal tower.

C.4.2 Other sources of error

The CFT matrix elements in the cost functions are computed using scaling dimensions and conformal spins extracted from energies and momenta of the excited states,

[TABLE]

where $A,B$ and $\gamma>1$ are non-universal constants. It then follows that the extracted scaling dimensions $\Delta_{\alpha}$ and the central charge $c$ have finite-size errors compared with $\Delta^{\textsl{\tiny CFT}}_{\alpha}$ and $c^{\textsl{\tiny CFT}}$. This error can be reduced by going to large sizes $N$ and through a large $N$ extrapolation. The scaling with $1/N$ depends on $\gamma$, which relies on specific irrelevant perturbations in the lattice Hamiltonian. For the Ising model $\gamma=3$, and we have reached several hundred spins in Zou et al. (2018) which makes this error on the order of $10^{-7}$ for $\Delta_{\epsilon}$ and $\Delta_{\sigma}$. This can be negligible compared to the truncation error in our range of system sizes $N\leq 48$.

There are also errors associated with numerical diagonalization of the Hamiltonian. Here we follow Zou et al. (2018) to use periodic uniform matrix product states (puMPS) as the diagonalization method, which was shown to result in energy eigensates with numerical errors that grow as we increase their energy. The fidelity of low energy eigenstates can be improved systematically by increasing the bond dimension of the puMPS. Here for the Ising model we use bond dimensions $D\leq 22$ for $N\leq 48$. Errors in fidelity of the first 23 eigenstates ($\Delta\leq 3+1/8$) are at most on the order of $10^{-6}$.

However, we note that the errors introduced during the diagonalization of the Hamiltonian will result in errors in $a_{\alpha}(N)$ that grow with both the scaling dimension $\Delta_{\alpha}$ of the operator and the system size $N$. This is because, in the cost function Eq. (77), the coefficient $\delta a_{\alpha}(N)$ is multiplied by a matrix element that scales as $N^{-\Delta_{\alpha}}$. To make this error always smaller than the truncation error, in this paper we only kept up to second level descendants in the set $\mathcal{A}$ and system sizes $N\leq 48$ when analysing the Ising model.

A more in-depth discussion of numerical errors in the Ising model can be found in the last appendix.

C.4.3 Error in numerical estimates of OPE coefficients

The OPE coefficients are approximately computed by

[TABLE]

where $\mathcal{O}_{\phi_{\beta}}$ is a lattice operator that corresponds to $\mathcal{O}^{\textsl{\tiny CFT}}_{\phi_{\beta}}\approx\phi^{\textsl{\tiny CFT}}_{\beta}$. Expanding in terms of scaling operators,

[TABLE]

where $a_{0}\approx 1$ and $\psi^{\textsl{\tiny CFT}}_{\beta^{\prime}}$ represents other scaling operators. Then

[TABLE]

We see that the error of $C_{\alpha_{\beta}\gamma}$ has two contributions. The first contribution, Eq. (83), contributes to a constant proportional to $a_{0}-1$, which is determined by the accuracy of the expansion coefficients of each lattice operator that are used to construct $\mathcal{O}_{\phi_{\beta}}$. The second contribution, Eq. (84), scales as $N^{-(\Delta_{\beta^{\prime}_{0}}-\Delta_{\beta})}$, where $\Delta_{\beta^{\prime}_{0}}$ is the scaling dimension of $\psi^{\textsl{\tiny CFT}}_{\beta^{\prime}_{0}}$ that appears in Eq. (84).

Therefore, in order to increase the accuracy of $C_{\alpha\beta\gamma}$, we can either compute Eq. (81) at larger sizes $N$, or obtain more significant digits of $a_{0}$. We shall see in the next appendix for the Ising model that both could lead to a significant improvement of accuracy.

Appendix D The critical Ising model

In this section we first exactly compute some matrix elements in the low energy spectrum of the Ising model using the free fermion representation. Then we use these exact matrix elements to obtain an exact expression for some of the (numerical) coefficients in Table 1. Finally, we analyse the numerical results.

D.1 Free fermion representation

Consider the critical Ising model with periodic boundary conditions,

[TABLE]

where site $j=N+1$ is identified with site $j=1$. The model has a $Z_{2}$ symmetry generated by

[TABLE]

It is easy to check that $[G_{S},H]=0$. $G_{S}$ and $H$ can be then simultaneously diagonalized, resulting in the eigenvectors of $H$ divided into parity even ($G_{S}=1$) and parity odd ($G_{S}=-1$) sectors.

The Jordan-Wigner transformation

[TABLE]

maps the Ising model with $N$ spins to a spinless fermion chain with $2N$ Majorana fermions, where

[TABLE]

and

[TABLE]

Local spin operators with odd $Z_{2}$ symmetry are mapped to a string of fermion operators, while those with even $Z_{2}$ symmetry are mapped to local operators in the fermion picture. We list some examples in table 2.

Let us represent the Ising Hamiltonian with the fermionic variables,

[TABLE]

One has to be careful with the boundary term. In the even $Z_{2}$ sector, the fermionic chain has the anti-periodic boundary condition, $\psi(2N+j)=-\psi(j)$, which is usually referred to as the Neveu-Schwarz (NS) sector. On the other hand, in the odd $Z_{2}$ sector, the fermionic chain has periodic boundary condition, $\psi(2N+j)=\psi(j)$, which is usually referred to as the Ramond (R) sector. We shall only consider the even $Z_{2}$ sector below.

Assuming the NS boundary condition, the Hamiltonian can be written more compactly as

[TABLE]

Note that the Hamiltonian is quadratic in fermionic variables. This makes it a free theory, which can be solved exactly using a Fourier transformation (below).

D.2 Symmetry and self-duality

In fermionic variables, the generator of the $Z_{2}$ symmetry can be expressed as

[TABLE]

It is then easy to see that $G_{S}$ commutes with fermionic bilinears, of the form $\psi(j)\psi(j^{\prime})$, but anti-commutes with operators linear in $\psi(j)$.

The Ising model at criticality possesses the famous Kramers-Wannier self-duality. This becomes a translation in the Majorana fermion picture,

[TABLE]

The Hamiltonian Eq. (94) is then manifestly invariant under the duality transformation. Note that applying the duality transformation twice corresponds to a translation by one site in spin variables.

Using the fermionic representation of local spin operators (Table 2), it is easy to see how they transform into each other under the duality transformation.

[TABLE]

We can then combine them into duality even operators (e.g., $XX+Z$) and duality odd operators (e.g., $XX-Z$). In the Ising model, operators that are $Z_{2}$ even and duality even belong to the conformal tower of the identity primary, while operators that are $Z_{2}$ even and duality odd belong to the $\epsilon$ tower. $Z_{2}$ odd operators (which are not considered here) belongs to the $\sigma$ tower.

D.3 Ground state correlation functions

Define the Fourier modes of fermion operators as

[TABLE]

where the boundary conditions impose

[TABLE]

The inverse Fourier transform is

[TABLE]

Eqs (90)-(91) imply

[TABLE]

It follows that $\psi(-p),\psi(p)$ can be understood as fermionic creation and annihilation operator of the mode $p$. Therefore we shall only include the modes with $p>0$ as independent variables.

The Hamiltonian can be rewritten as

[TABLE]

It follows that the ground state satisfies

[TABLE]

The two point correlation function in momentum space is then

[TABLE]

for $p_{n}>0$ and [math] otherwise. Fourier transforming back to position space yields

[TABLE]

In the thermodynamic limit $N\rightarrow\infty$, the above sum becomes an integral,

[TABLE]

and [math] if $j-l$ is even and nonzero.

For later use, we also present how the correlation function at finite $N$ behaves,

[TABLE]

if $j-l$ is odd.

When $j=l$, we have

[TABLE]

We then have all the two point correlation functions in the thermodynamic limit.

For example,

[TABLE]

And similarly,

[TABLE]

We can then derive higher point correlation functions by Wick’s theorem. For example, in the thermodynamic limit,

[TABLE]

Then

[TABLE]

The ground state expectation value of a lattice operator in the thermodynamic limit gives the coefficient of the identity operator in the corresponding CFT operator. For example

[TABLE]

where $\cdots$ represents other scaling operators. In this way, we obtain the coefficient in front of the identity operator for all $Z_{2}$ even operators in Table 1, as listed in Table 3.

D.4 Excited states

Excited states are created by applying creation operators $\psi^{\dagger}(p_{n})(p_{n}>0)$ on the ground state. There are two sets of creation operators at low energy, those with $p_{n}$ near $p=0$ and near $p=\pi$, corresponding to chiral and anti-chiral excitations. In the $Z_{2}$ even sector, there is an even number of fermions. The lowest lying excitations are

[TABLE]

The above phases will be determined shortly.

Matrix elements of lattice operators involving these excited states can then be computed by multi-point correlation functions of Majorana operator. For example,

[TABLE]

where the first equality follows from Table 2, the second equality follows from the Wick theorem, the third equality follows from Eq. (108) and Eq. (103), and the last two equalities follows from $p_{0}+p_{N-1}=\pi$.

At large sizes,

[TABLE]

As stated in Appendix C, we fix the phase $\theta_{\epsilon}=\pi/2$ by requiring the above matrix element to be real and positive.

Comparing to the CFT result

[TABLE]

and

[TABLE]

(derived from Eq. (45)), we can read off

[TABLE]

where $\cdots$ may contain $\partial^{2}_{x}\epsilon^{\textsl{\tiny CFT}}$ and higher scaling dimensions in the $\epsilon$ tower, as well as operators in the identity tower.

Similarly, we can compute

[TABLE]

where $s_{T}=2$ is the conformal spin of $T$.

Expanding it with respect to $1/N$ gives

[TABLE]

By requiring it to be negative, we fix $\theta_{T}=\pi$. Then we can compare it to

[TABLE]

where $c=1/2$ is the central charge, to obtain

[TABLE]

where $\cdots$ contains other scaling operators. In the same way,

[TABLE]

Now we can combine Eqs. (122),(130),(134),(135) to obtain

[TABLE]

where $\cdots$ contains other scaling operators with scaling dimension 3 or higher.

Since $XX$ and $Z$ are related by a duality transformation, they have the same coefficients in front of operators in the identity tower, and opposite coefficients in front of operators in the $\epsilon$ tower. Then

[TABLE]

Proceeding as above for other lattice operators as listed in Table 2, we reproduce part of Table 1 analytically, as listed in Table 3. We note that the subleading term in Eq. (112) is important in deriving the coefficient in front of $\partial^{2}_{\tau}\epsilon^{\textsl{\tiny CFT}}$ for $XIX$ and $ZZ$, which are quartic in fermionic variables. The coefficient in front of $\partial^{2}_{x}\epsilon^{\textsl{\tiny CFT}}$ cannot be computed by the matrix elements $\langle\epsilon|O(j)|0\rangle$ because $\langle\epsilon^{\textsl{\tiny CFT}}|\partial^{2}_{x}\epsilon^{\textsl{\tiny CFT}}|0^{\textsl{\tiny CFT}}\rangle=0$. Instead, we have to use matrix elements involving a state in the $\epsilon$ tower with non-vanishing conformal spin, such as $|\partial\epsilon\rangle$.

In Table 3, we also show the first $5$ digits of each analytically computed coefficient, to be compared with numerical results in the main text.

Comparing Table 3 with Table 1 in the main text, we see that in general, the coefficients in front of CFT operators with lower scaling dimensions are, as expected, more accurate.

We note in passing that, in order to reproduce the expansion for $Z_{2}$ odd operators analytically, we have to work with states in the Ramond sector and string operators in the fermion language. This is more complicated and we omit it here.

D.5 Analysis of numerical computations

Next we discuss the extrapolation to large system size $N$ of the Ising model of the numerical estimates for both the expansion coefficients and the OPE coefficients.

D.5.1 Convergence of expansion coefficients for the Ising model

We start with the numerical extrapolation of some of the expansion coefficients $a_{\alpha}$ presented in Table 1 in the main text. We show that the convergence of $a_{\alpha}(N)$ is as $1/N^{p}$, where $p=\Delta_{c}-\Delta_{\alpha}$ is the predicted power law in Eq. (78).

The first example is $\mathcal{O}=XZ+ZX$ with

[TABLE]

The error due to using a truncated set $\mathcal{A}$ of CFT scaling operators is determined by $\Delta_{c}>2+1/8$. It turns out that this particular operator does not have a contribution from the $\sigma$ tower at level $3$, and therefore $\Delta_{c}=4+1/8$. Therefore,

[TABLE]

The extrapolation of the numerical results is shown in Fig. 3.

The second example is $\mathcal{O}=ZZ$ and

[TABLE]

In this case, the leading complementary operator has scaling dimension $\Delta_{c}=4$ for the identity tower and $\Delta_{c}=5$ for the $\epsilon$ tower. Then

[TABLE]

The numerical results are shown in Fig. 4.

We have some additional comments on the extrapolation of the expansion coefficients for lattice operators. First, we cannot determine $\Delta_{c}$ a priori in general. Instead, we have to try extrapolation using different possible $\Delta_{c}$ to make the best fit.

Second, for some coefficient $a_{\alpha}(N)$, the error in numerical diagonalization may be important in the extrapolation. For example, this happens for $a_{\partial^{2}_{\tau}\sigma}\approx-0.017$ for $\mathcal{O}=X$ with Eq. (138), see Fig. 5.

Third, the extrapolation assumes the asymptotic scaling of $\delta a_{\alpha}(N)$ at large sizes. Numerically, we can determine $a_{\alpha}$ more accurately by using data from larger sizes, if other sources of error are negligible. In the operators that are considered, we find that for $\mathcal{O}=-i(YZ+ZY)$ with Eq. (138), the coefficients $a_{\partial^{2}_{\tau}\sigma}$ and $a_{\partial^{2}_{x}\sigma}$ are only obviously below $10^{-2}$ when extrapolating with data up to $N=96$, see Fig. 6.

D.5.2 OPE coefficients from the Ising model

According to Table 1 in the main text, we can three different ansatz for $\mathcal{O}_{\sigma}$, which all have their corresponding CFT operator in the form of Eq. (82), with $\phi^{\textsl{\tiny CFT}}_{\beta}=\sigma^{\textsl{\tiny CFT}},a_{0}\approx 1$. They are

[TABLE]

We quote the expansion coefficients that are used here for reader’s convenience,

[TABLE]

where we omit the $\partial^{2}_{x}\sigma^{\textsl{\tiny CFT}}$ term because it does not contribute to the OPE coefficient. In the following, we shall regard the coefficient of $\sigma^{\textsl{\tiny CFT}}$ in the above two expansions numerically the same, as they coincide with the highest accuracy (6 digits) among all coefficients that are computed.

In order to have $a_{0}\approx 1$, we determine

[TABLE]

Since $a_{0}$ has 6 significant digits and its error can be negligible to the finite-size errors below, we shall ignore the difference between $a_{0}$ and $1$.

The subleading operator in Eq. (82) for $\mathcal{O}_{\sigma 1}$ and $\mathcal{O}_{\sigma 2}$ is $\partial^{2}_{\tau}\sigma^{\textsl{\tiny CFT}}$, with coefficient $a_{(\partial^{2}_{\tau}\sigma)1}=-0.017\mu_{1}\approx-0.021$ and $a_{(\partial^{2}_{\tau}\sigma)2}=-0.820\mu_{2}\approx-1.02$. Therefore, the corresponding OPE coefficients are, according to Eqs. (83)-(84),

[TABLE]

where we have defined

[TABLE]

and where we have used Eq. (A.2). Substituting numbers into Eqs. (156)-(157) gives

[TABLE]

Eqs. (160)-(161) are confirmed with numerical results, see Fig. 7.

We see that using $\sigma^{\textsl{\tiny CFT}}\sim\mathcal{O}_{\sigma 1}$ results in much smaller errors, although they scale with the same power of $N$ as the errors in $\sigma^{\textsl{\tiny CFT}}\sim\mathcal{O}_{\sigma 2}$. This, of course, originates from the fact that the coefficient of $\partial^{2}_{\tau}\sigma^{\textsl{\tiny CFT}}$ in Eq. (153) has a much smaller amplitude than that in Eq. (154). The point of introducing the third lattice realization $\mathcal{O}_{\sigma 3}$ is to completely eliminate this contribution. Therefore, $\nu_{3}\approx-0.017/0.820\approx-0.021$ and it follows that $\mu_{3}\approx 1.27$. Ideally, this would result in

[TABLE]

However, since $\nu_{3}$ has not been determined with enough accuracy, the effect $\partial^{2}_{\tau}\sigma$ is not completely removed, and we did not observe a $1/N^{4}$ scaling. Instead, we find a $1/N^{2}$ scaling with an almost vanishing coefficient, see Fig. 8.

We can analyze the cases for the other OPE coefficient $C_{\sigma\epsilon\sigma}$ in the same way. We have used 4 different lattice realizations that approximately correspond to $\epsilon^{\textsl{\tiny CFT}}$,

[TABLE]

where $\mu_{1}\approx\mu_{2}=1/0.63662\approx 1.5708$ and $\mu_{3}=1/(-1.0876)\approx-0.92527$. The last operator is to eliminate the $\partial^{2}_{\tau}\epsilon^{\textsl{\tiny CFT}}$ contribution. Therefore, $\nu_{4}=-0.010/0.089\approx-0.11$ and $\mu_{4}\approx\mu_{1}/(1+\nu_{4})\approx 1.76$. However, since

[TABLE]

it does not improve the scaling of $C_{\sigma\epsilon\sigma}$. Therefore, we shall not show the OPE coefficient computed with $\mathcal{O}_{\epsilon 4}$, which almost coincides with that with $\mathcal{O}_{\epsilon 1}$. The numerical results are shown in Fig. 9

We see that only $\mathcal{O}_{\epsilon 3}$ has significant finite-size error in this OPE coefficient, which hints that in the expansion of $\mathcal{O}^{\textsl{\tiny CFT}}_{\epsilon 3}$ there is a significant contribution from $(L_{-4}+\bar{L}_{-4})\epsilon^{\textsl{\tiny CFT}}$. For the other two operators $\mathcal{O}_{\epsilon 1}$ and $\mathcal{O}_{\epsilon 2}$ , the above extrapolation suggests that the only subleading CFT operators are derivative descendants, which do not contribute to the OPE coefficient. This may be explained by the fact that they are quadratic in fermionic variables and the Ising model is a free fermion.