Spectrum of the Wilson-Fisher conformal field theory on the torus

TL;DR

This paper analyzes the finite-size spectrum of the Wilson-Fisher O(N) conformal field theory on a torus using epsilon expansion and numerical methods, revealing universal spectral features that identify different universality classes.

Contribution

It derives universal effective Hamiltonians for the zero modes and combines analytical epsilon expansion with numerical diagonalization to study the critical spectrum.

Findings

Effective Hamiltonians are N-dimensional anharmonic oscillators.

Numerical solutions match analytical predictions.

Spectral differences serve as fingerprints for universality classes.

Abstract

We study the finite-size spectrum of the O($N$) symmetric Wilson-Fisher conformal field theory (CFT) on the $d=2$ spatial-dimension torus using the expansion in $\epsilon=3-d$. This is done by deriving a set of universal effective Hamiltonians describing fluctuations of the zero momentum modes. The effective Hamiltonians take the form of $N$-dimensional quantum anharmonic oscillators, which are shown to be strongly coupled at the critical point for small $\epsilon$. The low-energy spectrum is solved numerically for $N = 1,2,3,4$. Using exact diagonalization (ED), we also numerically study explicit lattice models known to be in the O($2$) and O($3$) universality class, obtaining estimates of the low-lying critical spectrum. The analytic and numerical results show excellent agreement and the critical low energy torus spectra are qualitatively different among the studied CFTs, identifying them as a useful fingerprint for detecting the universality class of a quantum critical point.

Tables

Table 1. Table 1: Low-lying spectrum of the critical Ising model on the square torus from ϵ italic-ϵ \epsilon -expansion. The states shaded gray are odd under the global ℤ 2 subscript ℤ 2 \mathbb{Z}_{2} symmetry, while the unshaded states are ℤ 2 subscript ℤ 2 \mathbb{Z}_{2} even. Here we parametrize the momentum by κ = L ​ | k | / ( 2 ​ π ) 𝜅 𝐿 𝑘 2 𝜋 \kappa=L|k|/(2\pi) .

$\tau =i$	$\kappa =0$	$\kappa =1$	$\kappa =\sqrt{2}$	$\kappa =2$
$LE$	-0.12
	1.70
	5.04
		6.28
		8.76
			8.90
	8.90
			11.23
		12.29
				12.55

Table 2. Table 2: Low-lying spectrum of the critical O(2) model on the square torus from ε 𝜀 \varepsilon -expansion, including the ground state energy. These states are obtained from the effective Hamiltonian which gives the splitting of the k = 0 𝑘 0 k=0 Fock vacuum. The ℓ = 0 ℓ 0 \ell=0 states are non-degenerate while the ℓ > 0 ℓ 0 \ell>0 states are two-fold degenerate.

$\tau =i$	$\ell =0$	$\ell =1$	$\ell =2$	$\ell =3$	$\ell =4$
$LE$	-0.986
		0.96
			3.66
	4.74
				6.79
		8.31
					10.32
			12.12
				16.19
					20.52

Table 3. Table 3: Low-lying spectrum of the critical O(3) model on the square torus from ε 𝜀 \varepsilon -expansion, including the ground state energy. These states are obtained from the effective Hamiltonian which gives the splitting of the k = 0 𝑘 0 k=0 Fock vacuum. The states have degeneracy 2 ​ ℓ + 1 2 ℓ 1 2\ell+1 .

$\tau =i$	$\ell =0$	$\ell =1$	$\ell =2$	$\ell =3$	$\ell =4$
$LE$	-2.28
		-0.66
			1.74
	2.95
				4.59
		6.20
					7.81
			9.69
				13.45
					17.47

Table 4. Table 4: The coefficients of the strong-coupling expansion for the two-dimensional quartic oscillator.

$N=2$	$\ell$
	0	1	2	3	4
$c_{0,\ell ,0}$	1.47715	3.39815	5.65434	8.09067	10.7583
$c_{0,\ell ,1}$	0.258539	0.447039	0.605913	0.747439	0.877189
$c_{0,\ell ,2}$	-0.012345	-0.015633	-0.017109	-0.017109	-0.018469
$c_{0,\ell ,3}$	0.000903	0.000806	0.000697	0.000613	0.000548
$c_{1,\ell ,0}$	6.00339	8.70045	11.53475	14.50868	17.61616
$c_{1,\ell ,1}$	0.554312	0.682554	0.80824713	0.92837	1.04294
$c_{1,\ell ,2}$	-0.011291	-0.012457	-0.013714	-0.014705	-0.015489
$c_{1,\ell ,3}$	0.000133	0.000236	0.000293	0.000315	0.000320

Table 5. Table 5: The coefficients of the strong-coupling expansion for the three-dimensional quartic oscillator.

$N=3$	$\ell$
	0	1	2	3	4
$c_{0,\ell ,0}$	2.393644	4.478039	6.830308	9.401160	12.159017
$c_{0,\ell ,1}$	0.357801	0.529165	0.678421	0.813557	0.938665
$c_{0,\ell ,2}$	-0.014371	-0.016492	-0.017576	-0.018230	-0.018667
$c_{0,\ell ,3}$	0.000865	0.000749	0.000651	0.000578	0.000521
$c_{1,\ell ,0}$	7.335730	10.099944	13.004563	16.046193	19.217579
$c_{1,\ell ,1}$	0.618248	0.746036	0.869032	0.986315	1.09832
$c_{1,\ell ,2}$	-0.011790	-0.013117	-0.014238	-0.015120	-0.015816
$c_{1,\ell ,3}$	0.000188	0.000271	0.000307	0.000319	0.000320

Table 6. Table 6: The coefficients of the strong-coupling expansion for the four-dimensional quartic oscillator.

$N=4$	$\ell$
	0	1	2	3	4
$c_{0,\ell ,0}$	3.398150	5.624339	8.090668	10.758265	13.600878
$c_{0,\ell ,1}$	0.447038	0.605918	0.747451	0.877202	0.998248
$c_{0,\ell ,2}$	-0.015634	-0.017110	-0.017939	-0.018466	-0.018830
$c_{0,\ell ,3}$	0.000806	0.000697	0.000612	0.000547	0.000497
$c_{1,\ell ,0}$	8.700454	11.534729	14.508675	17.616152	20.849517
$c_{1,\ell ,1}$	0.682554	0.808247	0.928370	1.042942	1.152539
$c_{1,\ell ,2}$	-0.012457	-0.013714	-0.014704	-0.015489	-0.016111
$c_{1,\ell ,3}$	0.000236	0.000293	0.000315	0.000320	0.000319

Table 7. Table 7: Low-lying spectrum of the O(2) model from ED for geometry τ = i 𝜏 𝑖 \tau=i . The given values are obtained by averaging over the results for the different considered models/parameters. Only the lowest levels in the fully symmetric representation regarding the point-group symmetry are listed.

$\tau =i$	$\kappa =0$					$\kappa =1$
	$\ell =0$	$\ell =1$	$\ell =2$	$\ell =3$	$\ell =4$	$\ell =0$	$\ell =1$	$\ell =2$	$\ell =3$	$\ell =4$
	0						1.21
		0.24				1.39
			0.66					1.58
	1.00					1.79
$L\Delta /{\Delta }_{\ell =0}$				1.21			2.02
		1.62							2.12
					1.86					2.78
			2.40					2.92
				3.14
					3.79

Table 8. Table 8: Low-lying spectrum of the O(2) model from ED for geometry τ = 1 2 + 3 2 ​ i 𝜏 1 2 3 2 𝑖 \tau=\frac{1}{2}+\frac{\sqrt{3}}{2}i . The given values are obtained by averaging over the results for the different considered models/parameters. Only the lowest levels in the fully symmetric representation regarding the point-group symmetry are listed.

$\tau =\frac{1}{2}+\frac{\sqrt{3}}{2}i$	$\kappa =0$					$\kappa =1$
	$\ell =0$	$\ell =1$	$\ell =2$	$\ell =3$	$\ell =4$	$\ell =0$	$\ell =1$	$\ell =2$	$\ell =3$	$\ell =4$
	0						1.26
		0.25				1.46
			0.68					1.65
	1.00					1.83
$L\Delta /{\Delta }_{\ell =0}$				1.24			2.09
		1.64							2.19
					1.90					2.67
			2.39					2.85
				3.23
					3.84

Table 9. Table 9: Low-lying spectrum of the O(3) model from ED for geometry τ = i 𝜏 𝑖 \tau=i . The given values are obtained from the results of the O(3) bilayer model. Only the lowest levels in the fully symmetric representation regarding the point-group symmetry are listed.

$\tau =i$	$\ell =0$	$\ell =1$	$\ell =2$	$\ell =3$	$\ell =4$
	0
		0.29
			0.69
	1.00
$L\Delta /{\Delta }_{\ell =0}$				1.19
		1.55
					1.76
			2.17
				2.83
					3.28

Table 10. Table 10: Low-lying spectrum of the O(3) model from ED for geometry τ = 1 2 + 3 2 ​ i 𝜏 1 2 3 2 𝑖 \tau=\frac{1}{2}+\frac{\sqrt{3}}{2}i . The given values are obtained from the results of the O(3) bilayer model. Only the lowest levels in the fully symmetric representation regarding the point-group symmetry are listed.

$\tau =\frac{1}{2}+\frac{\sqrt{3}}{2}i$	$\ell =0$	$\ell =1$	$\ell =2$	$\ell =3$	$\ell =4$
	0
		0.24
			0.60
	1.00
$L\Delta /{\Delta }_{\ell =0}$				1.05
		1.49
					1.59
			2.08
				3.01
					3.58

Equations

$$H=\int d^{d}x\left(\frac{1}{2}{\Pi }_{\alpha }^2+\frac{1}{2}{\left(\nabla {\varphi }_{\alpha }\right)}^2+\frac{s_0}{2}{\varphi }_{\alpha }^2+\frac{u}{4!}{\varphi }_{\alpha }^4+\Lambda \right)$$

$$S_{d+1}=\frac{2}{\Gamma (d/2){\left(4\pi \right)}^{d/2}}$$

$$u=Z_4\frac{{\mu }^{ϵ}g}{S_{d+1}}$$

$$s_0=s_c+Z_{2}s$$

$$\Lambda =Z_{\Lambda }\frac{s^2{Z}_2^2{S}_{d+1}}{{\mu }^{ϵ}}$$

$${\varphi }_{\alpha }(x)$$

$${\Pi }_{\alpha }(x)$$

$$[{\phi }_{\alpha },{\pi }_{\beta }]=i{\delta }_{\alpha \beta },[b_{\alpha }(k),b_{\beta }^{†}(k^{\prime })]={\delta }_{\alpha \beta }{\delta }_{k{k}^{\prime }}$$

$$H$$

$$H_0$$

$${\mathcal{E}}_0={\mathcal{A}}^{d/2}\Lambda +\frac{N}{2}\sum _{k\ne 0}{\omega }_k,$$

$$V$$

$$H_0\Psi [{\phi }_{\alpha }]|k,\alpha ;k^{\prime },\beta ;\cdots ⟩=\left({\mathcal{E}}_0+{\omega }_k+{\omega }_{k^{\prime }}+\cdots \right)\Psi [{\phi }_{\alpha }]|k,\alpha ;k^{\prime },\beta ;\cdots ⟩$$

$$H_{\mathrm{eff}}|{\alpha }_0⟩=E_{\alpha }|{\alpha }_0⟩$$

$$E_{\alpha }={ϵ}_0+\mathcal{O}(V)$$

$$H_{\mathrm{eff}}={ϵ}_0{P}_0+P_{0}V{P}_0+P_{0}V\frac{1-P_0}{{ϵ}_0-H_0}V{P}_0+\mathcal{O}(V^3)$$

$$\sum _{a=1}^M{\Psi }_a[\phi ]{|\{k_i\}⟩}_a$$

$$\sum _{b=1}^M{H}_{\mathrm{eff},ab}{\Psi }_b[\phi ]|\{k_i\}{⟩}_b=E{\Psi }_a[\phi ]|\{k_i\}{⟩}_a$$

$$H_{\mathrm{eff}}$$

$$K$$

$$R$$

$$U$$

$$\phi \to {ϵ}^{-1/6}\phi ,\pi \to {ϵ}^{1/6}\pi$$

$$H_{\mathrm{eff}}={\mathcal{E}}_0+\frac{{ϵ}^{1/3}}{\sqrt{\mathcal{A}}}\left[K\frac{{\pi }^2}{2}+\frac{1}{2}R{ϵ}^{-2/3}{\phi }^2+\frac{U/ϵ}{4!}{\phi }^4\right]+\cdots$$

$$H_{\mathrm{eff}}={\mathcal{E}}_0+\frac{{ϵ}^{1/3}}{\sqrt{\mathcal{A}}}h(ϵ)$$

$$h(ϵ)$$

$$k$$

$$r$$

$$u$$

$$r=\frac{\mathcal{A}s}{{ϵ}^{2/3}}+{ϵ}^{1/3}\left(r_1+\mathcal{O}(ϵ)\right)$$

$$\mathcal{A}s=\mathcal{O}(ϵ)$$

$$V$$

$$\Psi [\phi ]|0⟩.$$

$${\mathcal{H}}_{\mathrm{eff},k=0}=h_{k=0}|0⟩⟨0|.$$

Full text

I Introduction

The identification of quantum critical behavior is an interesting problem in condensed matter and statistical mechanics. A major aspect of this is the emergence of universal low-energy behavior in the vicinity of a continuous quantum critical point, which is controlled by a CFT in the scaling limit.

In previous works, analysis of the torus spectrum has been used to study exotic deconfined quantum criticality, including transitions into phases with $\mathbb{Z}_{2}$ topological order Schuler et al. (2016); Whitsitt and Sachdev (2016) and conformal gauge theories coupled to fermionic matter Thomson and Sachdev (2016).

In this paper we explore the finite-size energy spectrum of the Wilson-Fisher CFT, also known as the critical O($N$) model, in $(d+1)$ spacetime dimensions. One case where the structure of the spectrum is well-understood is when the system is on the $d$-dimensional sphere $S^{d}$. In this case, conformal invariance implies the state-operator correspondence, which states that the energy spectrum takes the form $E_{n}=c\Delta_{n}/R$ where $c$ is the model-dependent speed of light, $R$ is the radius of the sphere, and the $\Delta_{n}$ are the scaling dimensions of the operators of the CFT in an infinite volume Cardy (1984, 1985). These scaling dimensions are extremely constrained by conformal invariance, and the operator spectrum of many interesting CFTs has been mapped out using methods such as exact results available in two spacetime dimensions Francesco et al. (1997) and the conformal bootstrap El-Showk et al. (2012, 2014); Kos et al. (2014a, b, 2016).

The state-operator correspondence has proven to be very useful in studying $(1+1)$-dimensional critical points, where numerically computing the spectrum on the circle is routinely done to accurately identify critical points Feiguin et al. (2007); Suwa and Todo (2015). However, in higher dimensions the curved geometry has proven to be difficult to implement accurately Alcaraz and Herrmann (1987); Weigel and Janke (2000); Deng and Blöte (2002); Brower et al. (2013, 2016). In light of these difficulties, it seems natural to instead study the universal energy spectrum on flat geometries such as the torus, where the energy spectrum still takes the form $E_{n}=c\xi_{n}/L$ for some universal set of constants $\xi_{n}$ dependent on the shape of the torus. However, the structure of the torus spectrum is not simply related to the operator content, so one must use perturbative field theory.

In this paper we will use the $\epsilon$-expansion, where $\epsilon=3-d$, to compute the critical energy spectrum. We note that this method has already been applied by us to compute the universal spectrum for the Ising CFT Schuler et al. (2016), but here we will generalize the result to general O($N$) models, and discuss how to compute the spectrum when the model is tuned away from criticality by a relevant perturbation. Two of us have also computed the spectrum on the torus directly in $(2+1)$-dimensions at leading order in the $1/N$ expansion Whitsitt and Sachdev (2016); the correspondence between the $\epsilon$- and large-$N$ expansions is discussed in Appendix E.

In this paper we will also numerically study the torus spectrum of explicit lattice models with critical points in the O(2) and O(3) universality classes using exact diagonalization (ED). We note that the numerical spectrum of the Ising CFT, O(1), has been discussed in detail in Ref. Schuler et al. (2016). We will show that the numerical computations show excellent agreement with the $\epsilon$-expansion results. Beyond that, we will demonstrate that the critical low-energy torus spectra are intrinsically different among the distinct CFTs considered in this paper, characterizing the interpretation of the low-energy critical torus spectrum as a universal fingerprint of the underlying CFT and as a useful tool for investigating quantum critical points with diverse methods.

In Section II we will introduce our model and method, discussing how to treat fluctuations of the zero mode non-perturbatively using an effective Hamiltonian method. The main technical section of the paper is Section III, which details the structure of the effective Hamiltonians and how to compute them using perturbative quantum field theory, giving examples for several important special cases to demonstrate how the method works in general. In Section IV we discuss how to numerically obtain the spectrum of the effective Hamiltonians in a few special cases, and give some results for the low-energy spectrum for the O($N$) models at $N=2,3$. Finally, in Section V we study several lattice models in the Wilson-Fisher universality class and numerically obtain their critical spectra using ED, and we compare numerical results to the analytic calculations.

II The method

The Wilson-Fisher CFT is described by the bare real-time Hamiltonian

[TABLE]

where the index $\alpha$ ranges from $1,...,N$. We are using the notation $\phi_{\alpha}^{2}\equiv\phi_{\alpha}\cdot\phi_{\alpha}$ and $\phi_{\alpha}^{4}\equiv\left(\phi_{\alpha}^{2}\right)^{2}$, so the model has full O($N$) symmetry. We suppress time-dependence, set the speed of light to unity, and note that fields satisfy the equal-time commutator $[\phi_{\alpha}(x),\Pi_{\beta}(x^{\prime})]=i\delta_{\alpha\beta}\delta^{d}(x-x^{\prime})$. We have included a bare ground state energy density $\Lambda$, which is needed to renormalize the ground state energy. The critical point is obtained by tuning $s_{0}=s_{c}$, while $u$ approaches a fixed value $u^{\ast}$. This is a strongly-coupled theory for any finite $d<3$ and $N$, but its universal properties can be computed as a power series in either $\epsilon=3-d$ or $1/N$.

We are interested in the finite-size spectrum of the above model on a spatial torus in $d=2$, which is parametrized by complex coordinates $x=x_{1}+ix_{2}$. We use the standard parametrization of the torus in terms of two complex periods, $\omega_{1}$ and $\omega_{2}$, and define the complex modular parameter $\tau\equiv\omega_{2}/\omega_{1}$ with real and imaginary parts denoted by $\tau=\tau_{1}+i\tau_{2}$. Below, we will often give results in terms of the length scale $L\equiv|\omega_{1}|$. The area of the torus is given by $\mathcal{A}=\mathrm{Im}\left(\omega_{2}\omega_{1}^{\ast}\right)=\tau_{2}L^{2}$. This geometry is pictured in Fig. 1.

In this paper we calculate the spectrum in the $\epsilon$-expansion, which forces us to introduce extra dimensions Wilson and Fisher (1972). To this end, we consider $d/2$ copies of the desired two-dimensional spatial geometry, which retains the point-group and modular symmetries of the system while avoiding the introduction of any additional unphysical parameters. We then expand in $\epsilon$, and set $\epsilon=1$ to obtain predictions for the $d=2$ system.

As usual, we will need to eliminate the bare couplings in favor of finite renormalized couplings. Because renormalization is entirely due to short-distance divergences, the finite-volume theory on a manifold with no curvature must have identical renormalization constants to the infinite volume theory, since the only new length scale is large compared to the cutoff. Here we use a modified minimal subtraction scheme Zinn-Justin (2002), where divergent poles in $\epsilon$ are subtracted with extra factors of

[TABLE]

attached. We introduce the renormalized coupling $g$ by

[TABLE]

In our calculations, we will always set $g$ to its fixed point value immediately after poles in $\epsilon$ have been subtracted. We also write

[TABLE]

where $s$ is a renormalized tuning parameter describing relevant perturbations across the critical point at $s=0$.

We are also interested in the dependence of the ground state energy as a function of the geometry of the torus and the tuning $s$. This requires introducing the counterterm

[TABLE]

which renders the energy density finite. In the above expressions, the $Z$ factors contain poles in $\epsilon$, and $\mu$ is an arbitrary energy scale. In principle one also needs to renormalize the fields $\phi_{\alpha}$, but these will not contribute to the leading order expressions so we ignore this here. For a review of these definitions and their relation to $L=\infty$ observables, see Appendix A.

The main technical feature of the $\epsilon$-expansion in a finite volume is the importance of the zero-momentum mode. Since the fields are gapless at the fixed point for $\epsilon\rightarrow 0$, the zero mode generates incurable infrared divergences in perturbation theory. These can be related to the failure of expanding around mean field theory, where the zero mode can have arbitrarily large fluctuations in the absence of interactions. This results in a continuum of zero mode excitations for the free field theory, whereas the interacting theory must have a gap and a discrete set of states, since the quartic term in the Hamiltonian will suppress fluctuations of the zero mode. Therefore, the free field theory with a zero mode is not the correct starting point in perturbation theory.

As was first realized by Lüscher Lüscher (1982), and further developed by others studying finite size effects in classical critical phenomena Brézin and Zinn-Justin (1985); Rudnick et al. (1985), the solution is to separate the zero mode and treat it non-perturbatively. Since the finite momentum modes have an effective gap, they can be safely integrated out in a path integral approach, leading to an effective action for the zero momentum modes which must be treated exactly. Here we pursue a Hamiltonian approach rather than a path integral approach, but our method is the same in principle. We note that a similar real-time approach was used to study low-energy spectrum of Yang-Mills theory on the torus in Ref. Lüscher (1983).

We expand the fields and their conjugate momenta as

[TABLE]

where $\omega_{k}=\sqrt{|k|^{2}+s_{0}}$ and $k\cdot x=\mathrm{Re}(kx^{\ast})$. The values of momentum summed over are determined by the shape of the torus, see Appendix C. Our expansion has been chosen such that the operators $\varphi_{\alpha}$, $\pi_{\alpha}$, and $b_{\alpha}(k)$ are dimensionless and have the commutation relations

[TABLE]

We now insert this expression into (1) and separate the Hamiltonian into a “free” and “interaction” part, where we insist that all zero-mode contributions are included in the interaction part.

[TABLE]

with the bare ground state energy

[TABLE]

and an interaction Hamiltonian,

[TABLE]

where $\chi_{\alpha}\equiv b_{\alpha}(k)+b_{\alpha}^{\dagger}(-k)$, and repeated latin indices are summed. We can now develop perturbation theory in $V$. The ground state energy (9) will be renormalized along with interactions in our final expressions and expanded to the appropriate order in $\epsilon$.

Since the zero mode does not appear in the unperturbed Hamiltonian, the unperturbed eigenstates are given by

[TABLE]

Here, the energies are determined by the Fock states created by the $b_{\alpha}^{\dagger}$ operators, but states can be multiplied by arbitrary normalizable functionals of the zero mode $\Psi[\varphi_{\alpha}]$. So the unperturbed states are infinitely degenerate, but this degeneracy is broken in perturbation theory.

We use a perturbation expansion developed by Bloch Bloch (1958) which is well-suited to dealing with degeneracy. The method involves deriving an effective Hamiltonian within each degenerate subspace whose spectrum gives the splitting of that subspace. So we find an operator $H_{\mathrm{eff}}$ such that

[TABLE]

where $|\alpha_{0}\rangle$ are the set of unperturbed degenerate states, and

[TABLE]

are the exact energies under the full interacting Hamiltonian, where $\epsilon_{0}$ is the unperturbed energy of the states $|\alpha_{0}\rangle$.

We review the derivation of Bloch’s effective Hamiltonian in Appendix B. The main result needed is that the effective Hamiltonian at leading order is given by

[TABLE]

where $P_{0}$ is the projection operator onto the degenerate subspace $|\alpha_{0}\rangle$. The calculation of the effective Hamiltonian will result in UV divergences due to summations over infinitely large momenta implicitly contained in Eq. (14), so this is the step where the theory is renormalized.

From Eq. (11), we can infer the action of $H_{\mathrm{eff}}$ on our degenerate subspaces. An arbitrary state in a given degenerate manifold takes the form

[TABLE]

where $M$ is the number of Fock states with the same energy. This degeneracy between inequivalent Fock states will be due to O($N$) symmetry or discrete rotation symmetry. The effective Hamiltonian will take the form of an $M\times M$ matrix which acts as

[TABLE]

This effective Hamiltonian will be a function of $\pi$ and $\varphi$. Due to the commutation relations, the action of $\pi$ on $\Psi_{a}[\varphi]$ is $\pi_{\alpha}=-i\frac{\partial}{\partial\varphi_{\alpha}}$, so this will be a set of $M$ coupled differential equations. In practice, symmetries of the interaction will allow us to consider smaller block-diagonal subspaces separately. For our model, the interaction conserves momentum and O($N$) rotations, which further constrains the form of $H_{\mathrm{eff}}$ and its eigenvectors.

III Calculation of the effective Hamiltonians

III.1 Structure of the effective Hamiltonians

Before proceeding with explicit calculations, we first describe the general structure of the effective Hamiltonian and its dependence on $\epsilon$, and discuss the perturbative spectrum of the Hamiltonian for small $\epsilon$. We give this discussion prior to explicit calculations because we will see that the perturbation theory is reordered in the scaling limit, leading to a modified expansion in fractional powers of $\epsilon$. We will find that the terms in $V$ do not contribute to the $\epsilon$-dependence of the spectrum that would be naïvely inferred by Eq. (10), so our analysis will aid us in correctly finding the leading contributions when we turn to explicit calculations. We will also highlight how the behavior of the spectrum changes depending on the magnitude of the tuning parameter $s$.

Here we will consider the effective Hamiltonian when the Fock state is non-degenerate ($M=1$ in Eq. (16)), since the analysis is similar in the general case. From Eq. (14), the effective Hamiltonian at leading order takes the form

[TABLE]

where the ellipses will include additional operators which will not appear until $\mathcal{O}(\epsilon^{2})$. The coefficient of each term will be a regular series in $\epsilon$. Here we write $r_{1}$ as the $\mathcal{O}(\epsilon)$ coefficient of the operator $\varphi_{\alpha}^{2}$. We have used the relations in Appendix A to set the quartic coupling to its fixed point.

Since we are interested in the critical regime, we first consider the theory for $s=0$. In this case, the coefficients of the quadratic and quartic terms are both $\mathcal{O}(\epsilon)$, but the structure of the spectrum can be made more clear by making the canonical transformation

[TABLE]

after which the Hamiltonian is given by

[TABLE]

We see that when $s=0$, the Hamiltonian takes the form

[TABLE]

where

[TABLE]

We see that $h(0)$ is a pure quartic anharmonic oscillator whose spectrum gives the spectrum of $H_{\mathrm{eff}}$ to order $\epsilon^{1/3}$. Furthermore, the leading corrections are given by obtaining the spectrum of $h(\epsilon)$ in a perturbation series in $\epsilon^{1/3}$. Thus, the spectrum of the Wilson-Fisher fixed point at small $\epsilon$ has mapped to the spectrum of a quartic anharmonic oscillator in a strong-coupling expansion.

If we repeat the above analysis for $s\neq 0$, we find that the quadratic coefficient of the reduced Hamiltonian $h(\epsilon)$ is modified to

[TABLE]

For $\mathcal{A}s\gtrsim\epsilon^{2/3}$, our previous analysis no longer holds. The spectrum will be given by a weak-coupling expansion around a simple harmonic oscillator Hamiltonian provided $\mathcal{A}s\gg\epsilon$. This is sensible because our entire approach has been based on the gaplessness of the zero mode, whereas a nonzero $s$ contributes a gap. For large enough values of $\mathcal{A}s$, we do not need to separate the zero mode, and we could have done a normal expansion around Gaussian field theory, which is equivalent to the weak coupling expansion in the current effective Hamiltonian approach. The crossover from a weakly-coupled oscillator with a particle-like spectrum to a strongly-coupled oscillator signals the breakdown of particle-like excitations at the quantum critical point. It is interesting that this occurs already for arbitrarily small values of $\epsilon$, reflecting the importance of interactions in confining the zero mode.

In this paper it will be desirable to use the same strong-coupling expansion of $h(\varepsilon)$ at $s=0$ and $s\neq 0$, so we will always assume

[TABLE]

or smaller.

The reordering of the perturbative expansion requires a modified analysis of our perturbation theory. This can be most easily visualized by representing terms in the effective Hamiltonian diagrammatically, associating extra factors of $\epsilon^{-1/6}$ with factors of $\varphi$ anticipating the utility of the transformation in Eq. (18). With this in mind, we can schematically rewrite the interaction Hamiltonian Eq. (10) after the canonical transformation to identify the individual terms with the correct powers of $\epsilon$:

[TABLE]

Here, the $M^{(n)}$ coefficients involve $n$ factors of finite-momentum modes. With this form of the interaction, we write down the vertices associated with $V$ in Fig. 2. Then in calculating the effective Hamiltonian from Eq. (14) we organize the $\epsilon$-expansion diagrammatically using these vertices. In practice it is easier to work directly with Eq. (10) to compute the effective Hamiltonian, but the correct order of each term’s contribution to the energy spectrum will be given by the $\epsilon$ coefficient pictured in Fig. 2.

III.2 Effective Hamiltonians for low-lying states

In this section we will give the explicit derivation of the effective Hamiltonians for the lowest-lying states in the Fock spectrum. We will perform the calculation for increasingly complex cases, with each example having an added subtlety compared to the previous case, after which the general structure for the effective Hamiltonian splitting an arbitrary Fock state should follow.

III.2.1 Fock vacuum

We begin by considering the splitting of the Fock vacuum. This will give us the lowest-lying zero-momentum states, including the energy gap. The unperturbed eigenstate is

[TABLE]

Since $P_{0}=|0\rangle\langle 0|$, the effective Hamiltonian will be of the form

[TABLE]

From Eq. (16), the Schrödinger equation acting on the unperturbed subspace reduces to

[TABLE]

where, using Eq. (14), the reduced Hamiltonian $h_{k=0}$ at one-loop is given by

[TABLE]

At this point we note that every term appearing in $h_{k=0}$ can be associated with a diagram. The three terms in this equation correspond to diagrams with zero, one, and two vertices respectively. Because the interaction $V$ conserves momentum, each vertex must also enforce momentum conservation. The expectation values and sums over $k$ implies one must contract all $\chi(k)$ propagators. Finally, the presence of the projector in the last term means we must contract the two vertices, preventing any disconnected diagrams from appearing.

Writing the effective Hamiltonian as

[TABLE]

we collect the one-loop terms which contribute to the effective Hamiltonian in Fig. 3. At $\mathcal{O}(\epsilon^{5/3})$ we encounter a nontrivial two-loop diagram, also pictured in Fig. 3, so we truncate the spectrum to order $\epsilon^{4/3}$.

Using Eq. $(\ref{intH})$, we obtain

[TABLE]

At this point, we need to evaluate these infinite sums in $3-\epsilon$ dimensions, and we need to renormalize the theory. These technical details are treated at length in the Appendices. In Appendix C it is shown how to evaluate infinite sums in arbitrary dimension, and in Appendix D these three expressions are explicitly evaluated, and the cancellation of all divergences is demonstrated. We obtain the following effective Hamiltonian:

[TABLE]

Thus, the lowest states in the spectrum of the O($N$) Wilson-Fisher fixed point are given by solving the quantum mechanics problem of an isotropic, $N$-dimensional anharmonic oscillator. Here, the special functions $f^{(d)}_{\nu}(\tau,s,\mu)$ are given explicitly in Eqns. (123-125), and the function $f^{(3-\epsilon)}_{-1/2}(\tau,s,\mu)$ should be expanded to first order in $\epsilon$.

The functions $f^{(d)}_{-1/2}$ and $f^{(d)}_{1/2}$ depend on the renormalization scale $\mu$, which can be eliminated by applying renormalization conditions on the $s\neq 0$ ground state energy and the energy gap. In this paper we will not eliminate $\mu$, since the infinite volume quantities are non-analytic around the critical point. However, our assumption that $s\sim\mathcal{O}(\epsilon)$ allows us to set $s$ to zero in many of the terms. We furthermore note that at the critical point, $s=0$, the $\mu$-dependence drops out and the spectrum is a universal function of $\tau$, $N$, and $\epsilon$. The coefficients of the Hamiltonian are modular invariant at $s=0$, which follows from Eq. (116) and the modular invariance of $\mathcal{A}=\tau_{2}L^{2}$.

After performing a canonical transformation similar to Eq. (18), this Hamiltonian takes the form

[TABLE]

We will primarily work with this form of the Hamiltonian in Section IV.

III.2.2 Single particle Fock states

We now consider the splitting of the single particle state

[TABLE]

where we assume there are no multi-particle Fock states with the same momentum and energy, so we only need to consider an $N$-fold degenerate manifold. This assumption should hold for the smallest values of $|k|$. The effective Hamiltonian can be written as a matrix equation acting on the vector of functions $\Psi_{\alpha}$:

[TABLE]

where the effective Hamiltonian $h_{k,\alpha\beta}$ can be represented by an $N\times N$ matrix whose components are

[TABLE]

and $|k,\alpha\rangle=b^{\dagger}_{\alpha}|0\rangle$.

The terms in this equation can also be given a diagrammatic representation: one can consider all diagrams with two external finite-momentum lines carrying momentum $k$. However, the external momenta do not need to be contracted with the vertices in any of the terms. As a consequence, we obtain a term

[TABLE]

simply by taking the diagrams in Fig. 3 and drawing a disconnected solid line in them. In addition to these, we will get new contributions due to the extra diagrams pictured in Figure 4. Only one of these diagrams contains a loop, so there will only be a single new divergence which is cancelled by the mass renormalization of $s_{0}$ in the first term of Eq. (35).

An explicit calculation similar to the one done for the Fock vacuum results in

[TABLE]

For $N>1$, this is an anisotropic anharmonic oscillator.

III.2.3 Two-particle Fock states

We now consider mixing between states of the form

[TABLE]

We have $|k_{1},\alpha;k_{2},\beta\rangle=S_{12}b^{\dagger}_{\alpha}(k_{1})b^{\dagger}_{\beta}(k_{2})|0\rangle$ where we require a Bose symmetry factor: $S_{12}=\frac{1}{\sqrt{2}}$ if $k_{1}=k_{2}$ and $\alpha=\beta$, and $S_{12}=1$ otherwise. We will assume that there are no one-particle or $\geq 3$-particle states with overlap with this state. However, we will consider the case where there exists inequivalent states $|k_{1},\alpha;k_{2},\beta\rangle$, $|k_{3},\alpha;k_{4},\beta\rangle$, such that

[TABLE]

Because the interaction $V$ conserves momentum, these states will only mix if

[TABLE]

Such states can contribute to the low-energy spectrum on the torus. For example, on the square ($\tau=i$) torus the states $|2\pi/L,\alpha;-2\pi/L,\beta\rangle$ and $|2\pi i/L,\alpha;-2\pi i/L,\beta\rangle$ are inequivalent but can mix.

The effective Hamiltonian is now calculated in a similar manner to the previous two cases, and there is an obvious diagrammatic generalization of the previous rules. We now draw diagrams with four external lines with momenta $k_{i}$, $i=1,...,4$. We then consider all possible contractions with either zero, one, or two vertices. We once again find a piece proportional to $h_{k=0}$, which involves the diagrams in Fig. 3 but with the four external lines contracted and disconnected to the vertices. In addition, we get contributions which are simply the diagrams in Fig. 4 but with a single additional finite-momentum line disconnected from the rest. Finally, we obtain the additional diagram pictured in Fig. 5, which can connect the inequivalent states considered above. This does not contain a loop, so it is finite.

An explicit calculation gives

[TABLE]

The first two lines of Eq. (41) can be given in terms of the the zero-particle and single-particle Hamiltonians, while the last term is new and contributes a constant shift in the energy. It is the last term which can mix two states unrelated by O($N$) symmetry, leading to a multidimensional Hamiltonian even when $N=1$.

III.2.4 Mixing between one- and two-particle Fock states

We now consider the effective Hamiltonian which couples the states $|k,\alpha\rangle$ and $|k_{1},\gamma;k_{2},\delta\rangle$. In order for these states to mix, we need

[TABLE]

While this set of equations is very restrictive, when $s=0$ it can be satisfied by choosing $k_{1}$ and $k_{2}$ to be collinear, so in our expressions where $\mathcal{A}s\lesssim\epsilon$, we always have $s\ll|k|^{2}$ and these states will mix. We will assume here that the only other mixing is due to O($N$) symmetry.

We can write the Hamiltonian mixing these two states as

[TABLE]

where the diagonal Hamiltonians are given by the one-particle and two-particle cases above, and the off-diagonal elements are

[TABLE]

where $P_{k}$ is the projector onto the degenerate subspace. We want all diagrams with three external lines with momenta $k$, $k_{1}$, and $k_{3}$. This leads to the diagram in Figure 6, which is finite, giving

[TABLE]

This coupling between Fock states with different particle number is another manifestation of the breakdown of well-defined quasiparticles at a quantum critical point.

At this point we hope it is clear how to generalize these results to more complicated cases. As one moves higher in the spectrum, more degeneracies between Fock states are possible and more complicated effective Hamiltonians will be needed.

IV Numerical solution of the effective Hamiltonians

It remains to find the spectrum of the effective Hamiltonians perturbatively in $\epsilon$. This must be done numerically, since the $\epsilon$-expansion is equivalent to a strong-coupling expansion of a quartic anharmonic oscillator, which has no analytic solution. Given that there are many effective Hamiltonians, and each needs to be solved to at least third order in perturbation theory to obtain the $\epsilon^{4/3}$ contribution, this becomes one of the biggest barriers to obtaining the spectrum.

Here we show how to obtain the spectrum for most of the low-lying effective Hamiltonians in the $N=1$ case, where the Hamiltonians are often one-dimensional quartic oscillators. For $N>1$, we will confine ourselves to solving the effective Hamiltonian which splits the Fock vacuum, where rotational invariance allows us to map it to a one-dimensional differential equation in “radial” coordinates.

In Appendix E we discuss of the $N=\infty$ limit, where we can solve the spectrum of the Fock vacuum effective Hamiltonian exactly. We show agreement with the results of the large-$N$ expansion in Ref. Whitsitt and Sachdev (2016) wherever possible.

IV.1 $N=1$

The results of this section were used in a previous publication by us Schuler et al. (2016). We exploit the fact that most of the effective Hamiltonians for the low-lying states can be written in the form (see e.g. Eq. (32))

[TABLE]

where

[TABLE]

Here, $\mathcal{E}$ is a known function of $\epsilon$, and

[TABLE]

for some constants $r_{1}$, $u_{1,2}$ which depend on the specific effective Hamiltonian in question. Then we need the spectrum of $h[r]$ in a power series in $r$. After writing

[TABLE]

expansion of the spectrum in $\epsilon$ can be written

[TABLE]

Assuming we know the coefficients $c_{n,m}$ which appear in this expansion, we can simply write down the spectrum to the desired order:

[TABLE]

Fortunately, the numerical calculation of the coefficients $c_{n,m}$ has been extensively studied in the literature. In Ref. Skála et al. (1999), Tables 7 and 10, the coefficients $c_{n,m}$ for $m=1,..,10$, $n=1,...,10$ are given with at least five digit accuracy. Thus, for one-dimensional effective Hamiltonians, the spectrum can be obtained.

However, there can still be multi-dimensional effective Hamiltonians for $N=1$ due to mixing between different Fock states. Fortunately, the form of the mixing is often very simple. In computing the low-lying spectrum, the first mixing one comes across is between two inequivalent two-particle states, as described in Section III.2.3. There, the effective Hamiltonian acting on these two states takes the form

[TABLE]

where $h_{2k,2}$ is of order $\epsilon$ and is independent of the zero mode. In this case, we can diagonalize the Hamiltonian by inspection, obtaining the energy splitting

[TABLE]

where $E_{2k,1}$ are the energies of $h_{2k,1}$.

Finally, we do encounter the states described in Section III.2.4 which mix the single-particle and two-particle Fock states. In this case, we must numerically calculate the contribution. However, this is made easier by the fact that the order $\epsilon^{5/6}$ term arises in first-order perturbation theory, so we only need the zeroth order wave functions. That is, we first numerically calculate the zeroth order wave functions which diagonalize $h_{k}$ and $h_{2k}$, and then we compute the overlap

[TABLE]

We calculate the unperturbed wave function numerically using the same shooting method described in Appendix F for the $N>1$ case.

The low-lying spectrum of the critical Ising model is given in Table 1.

IV.2 $N>1$

Here we will focus on the splitting of the Fock vacuum, where the Hamiltonian is an isotropic $N$-dimensional oscillator. We begin with Eq. (32):

[TABLE]

Then defining

[TABLE]

and the reduced couplings

[TABLE]

we obtain Eq. (51) except the coefficients $c_{n,m}$ will depend on $N$.

To make further progress, we take advantage of the fact that $h_{N}[r]$ can be written in spherical coordinates, after which it separates into a known angular equation and a one-dimensional radial equation. We go to hyperspherical coordinates

[TABLE]

Then the Laplacian can be written in the separable form

[TABLE]

where $\nabla_{S^{N-1}}^{2}$ is the Laplacian on the $(N-1)$-sphere. We will not use the coordinate representation of this operator, but instead use what is known of its spectrum Stein and Weiss (1971). The eigenvectors and eigenvalues are given by

[TABLE]

where the eigenfunctions can be given in terms of Gegenbauer polynomials of the $\cos\theta_{i}$, and the indices can range from

[TABLE]

The spectrum does not depend on the the $\ell_{i}$. This gives a degeneracy for a given eigenvalue $\ell$ of

[TABLE]

for $N\geq 3$. For $N=2$, there is only one eigenfunction for each $\ell$, but the states $Y_{\ell}$ and $Y_{-\ell}$ are degenerate (these are simply the states $e^{i\ell\theta}$ and $e^{-i\ell\theta}$), so the degeneracy is

[TABLE]

Finally, we note that the eigenfunctions are in the symmetric traceless tensor representation of O($N$), and we can label these representations using the eigenvalue $\ell$.

With these eigenfunctions, we can express our functionals as

[TABLE]

and the eigenvalue equation becomes

[TABLE]

Here we have introduced a radial quantum number $n$, which corresponds to the number of zeros in $R$. We have reduced our problem to a one-dimensional eigenvalue equation, and we wish to find the spectrum perturbatively in $r$. In analogy with the $N=1$ case, we write this expansion as

[TABLE]

We have obtained the coefficients of the perturbative expansion in $r$ for $N=2,3,4$, $\ell=0,...,4$ and $n=0,1$ numerically. We obtained these by first solving the $r=0$ equation numerically for the wave function and energy. We then used logarithmic perturbation theory Au and Aharonov (1979); Aharonov and Au (1979), which is well-suited to this problem because it allows one to find the coefficients of the expansion directly from the unperturbed energy and wave function without needing excited states. We discuss our approach in Appendix F, and give the coefficients of the expansion in Tables 4-6. Using these results, the energy to leading order is

[TABLE]

We give the lowest-lying states for $N=2$ and $N=3$ at the critical point in Tables 2 and 3.

V Numerical calculation of critical torus spectra from lattice models

In this section we investigate the critical torus spectra numerically using Exact Diagonalization (ED). We study explicit lattice models known to exhibit a transition in the O($N$) universality class, compute their spectrum on finite clusters and extrapolate it to the thermodynamic limit. We then compare these numerical results to the analytic calculations from the $\epsilon$-expansion. In the following, we will present results for $N=2$ and $N=3$ and refer to Ref. Schuler et al. (2016) for a detailed discussion about the $N=1$ (Ising) CFT.

V.1 $N=2$

We first consider the O(2) critical point, also known as the quantum XY model. Two different lattice models will be utilized to study the critical spectrum of this universality class numerically. The first is a spin-1 model, $S=1$, with single-ion anisotropy Zhang et al. (2013); Stoudenmire et al. (2014)

[TABLE]

We set the energy scale by choosing $J=1$. For small $D$ the system orders ferromagnetically in the $x$-$y$ spin plane, while for large $D$ the system approaches a product state of single spins with $S^{z}=0$. The phase transition is found to be in the XY/O(2) universality class. The parameter $J_{z}$ can be tuned within a range around zero to check the stability of our results.

The second model we consider is the spin-1/2 XY-bilayer model Stoudenmire et al. (2014); Helmes and Wessel (2015). It consists of two usual ferromagnetic XY layers with additional XY couplings between them. We denote a spin located on site $i$ in the first (second) layer as $S_{i,1(2)}$, the model is then described by

[TABLE]

We only consider positive couplings $J,J_{\perp}>0$ here and set our energy scale $J=1$. For large $J_{\perp}$ the system is described by a product state of singlets on each interlayer bond, whereas a XY-ferromagnet is formed in each plane for small $J_{\perp}$. The two phases are separated by a XY quantum critical point at $J_{\perp}=J_{\perp}^{c}=5.460(1)J$ Helmes and Wessel (2015)111The energy spectrum of this model is conserved under changing the sign of $J\rightarrow-J$ when the momenta for odd $S^{z}$ levels are shifted as $\mathbf{k}\rightarrow\mathbf{k}+(\pi,\pi)$. Therefore, the critical coupling $J_{\perp}^{c}$ is identical for both, ferromagnetic and antiferromagnetic $J=\pm 1$.

To calculate the critical torus energy spectrum, we compute the spectrum of the Hamiltonians at criticality on finite-size toric clusters numerically. We multiply the finite-size spectra with the linear system size $L$ to get rid of the dominant scaling and then extrapolate these spectra in $1/N$ to the thermodynamic limit. Further details about the numerical approach can be found in App. G.

In Fig. 7 we present the numerically obtained critical O(2) torus energy spectrum in the $\kappa=0$ and $\kappa=1$ sectors. Here and in the following we only show the lowest energy levels for $\ell\leq 4$ which are in the fully symmetric representation regarding the lattice point-group symmetry. We also restrict our discussion to square lattices, $\tau=i$. Results for triangular geometry $\tau=\frac{1}{2}+\frac{\sqrt{3}}{2}i$ are given in App. G. The spectrum is normalized such that the lowest gap in the $\ell=0$ sector is set to $\Delta_{\varepsilon_{T}}\equiv\Delta_{\ell=0}=1$. To demonstrate the stability of the numerical results and the universality of the spectrum, the different models and parameters considered for O(2) universality are shown with different symbols and colors in the plot. For $\kappa=0$ we also plot the $\epsilon$-expansion results with empty diamonds to compare them with the numerics. We find that they agree reasonably well and show a qualitatively identical structure. This highlights, that the critical torus energy spectrum is a universal fingerprint of the underlying CFT which is available from a wide variety of analytical and numerical approaches.

Here we also want to note, that the four relevant fields in the O(2) CFT Pelissetto and Vicari (2002) correspond to the lowest $\ell=1$, $\ell=2$ and $\ell=3$ levels as well as the first excited $\ell=0$ level in the critical torus spectrum (all $\kappa=0$). Interestingly, these are the four lowest states in the spectrum. Also in the Ising CFT, the two relevant fields correspond to the two lowest levels in the critical spectrum Schuler et al. (2016). This indicates, that it might be a general feature that relevant fields of the CFT have light analogues in the critical torus spectrum.

The numerical values for the critical torus spectrum are tabulated in App. G in Tab. 7 for the square geometry $\tau=i$ and in Tab. 8 for the triangular geometry $\tau=\frac{1}{2}+\frac{\sqrt{3}}{2}i$.

V.2 $N=3$

In this section, we study the critical torus spectrum of the O(3) CFT numerically. To do so, we again consider two different lattice models with a critical point known to be described by the O(3) universality class. The first model is the prototypical Heisenberg bilayer model Sandvik and Scalapino (1994); Wang et al. (2006); Lohöfer et al. (2015); Nishiyama (2016). It consists of spin-1/2 on two layers with nearest-neighbour Heisenberg intraplane couplings and Heisenberg interplane couplings on the rungs:

[TABLE]

$S_{i,n}$ denotes a spin on site $i$ in layer $n$. We set ferromagnetic intraplane couplings $J=-1$ and antiferromagnetic rung couplings $J_{2}>0$. For large $J_{2}$ the groundstate is a product state of singlets on each rung, whereas the groundstate for small $J_{2}$ is the direct product of Heisenberg ferromagnets within each plane. The phases are separated by an O(3) critical point.

The second model we want to investigate here is a Heisenberg model on a 2D square lattice with columnar dimerization of bonds Wenzel and Janke (2009); Matsumoto et al. (2001). The Hamiltonian for this ladder model is

[TABLE]

Every second horizontal bond on the lattice is chosen to be in the family $\langle i,j\rangle^{\prime}$, such that these dimerized bonds form ladders and every spin is part of exactly one dimerized bond. We set all couplings antiferromagnetic and set the energy scale $J=1$. For $J_{2}/J=1$ a Néel AFM is stabilized on a square lattice and for large $J_{2}/J$ a product state of singlets on the bonds $\langle i,j\rangle^{\prime}$ is formed as a groundstate. These phases are separated by an O(3) transition at the critical coupling $(J_{2}/J)_{c}=1.9096(2)$ Wenzel and Janke (2009).

We proceed similar to the case of $N=2$ to compute the critical torus energy spectrum for the O(3) CFT numerically. The critical spectrum is shown in Fig. 8 together with the $\epsilon$-expansion results. We again observe a decent agreement between the two methods and qualitatively identical critical torus spectra. Larger deviations for the second level in $\ell=4$ are probably related to difficulties in the extrapolation to the thermodynamic limit, as the available system sizes are strongly limited for these models. Although the critical torus spectra in the $\kappa=0$ sector seems to look very similar for the O(2) and O(3) CFTs, their degeneracy structure is inherently different. For O(3) the levels are $(2\ell+1)$-fold degenerate, whereas they are 2-fold (1-fold) degenerate for $\ell>0$ ($\ell=0$) in the O(2) CFT. The degeneracy structure is thus an important feature, yielding qualitatively very distinct critical torus spectra.

The numerical values for the critical torus spectrum of the O(3) CFT from ED are tabulated in App. G in Tab. 9 for the square geometry $\tau=i$. In Tab. 10 we additionally list the energy levels for the triangular geometry $\tau=\frac{1}{2}+\frac{\sqrt{3}}{2}i$.

VI Conclusions

In this paper, we have demonstrated how to compute the torus energy spectrum of the Wilson-Fisher CFT from quantum field theory in the $\epsilon$-expansion, showing the expected emergence of a spectrum determined entirely by the universality class of the CFT and the spatial geometry of the torus. We have also computed the spectrum of several explicit lattice models at their respective quantum critical points using ED, and have shown that the analytic and numerical calculations agree well with each other, highlighting that the finite-size spectrum is a useful universal fingerprint for identifying quantum criticality. In Fig. 9 we show a comparison of the critical torus low-energy spectra for $N=1,2,3$ in the $\kappa=0$ sector which substantiates this interesting aspect of a universal fingerprint. Additionally, we compare the critical torus spectra with operator scaling dimensions of the corresponding CFTs from literature El-Showk et al. (2012); Hasenbusch and Vicari (2011). We note a quantitative match between the low-energy critical torus spectrum at $\kappa=0$ and the (rescaled) operator dimensions. This match may be a coincidence, as there is no known mapping between the torus spectrum and the operator dimension. It is not known what constitutes the complete data for a three-dimensional CFT, and whether the torus spectrum is new data or whether is can be related to the set of operator scaling dimensions.

We hope that our work has shed light on the nature of the finite-size spectrum in CFTs where conformal invariance does not lead to simple and exact results, and that our calculations aid in identifying critical behavior in numerical studies of quantum lattice models by investigation of the critical energy spectrum. The methods used in this work should also be applicable to computing the finite-size spectra of other CFTs, such as an $\epsilon$-expansion approach to the Gross-Neveu-Yukawa CFT.

Another application of the present results is to a class of deconfined quantum critical points known as the O($N)^{\ast}$ models. These are described by the same field theory as the O($N$) model except that the order parameter carries a $\mathbb{Z}_{2}$ gauge redundancy, imposing an equivalence relation $\phi\sim-\phi$ Jalabert and Sachdev (1991); Chubukov et al. (1994). In these models, the order parameter $\phi$ can be thought of as the spinon in a $\mathbb{Z}_{2}$ spin liquid, whose condensation results in a confinement phase transition accompanied by broken symmetry associated with the internal O($N$) quantum number.

The prescription for relating our current results to the O($N)^{\ast}$ model was detailed in Refs. Schuler et al. (2016); Whitsitt and Sachdev (2016), and amounts to a selection rule on the allowed $\mathbb{Z}_{2}$ parity of the spectrum and on the inclusion of states with twisted boundary conditions $\phi(x+iy+n_{1}\omega_{1}+n_{2}\omega_{2})=\pm\phi(x+iy)$. The spectrum of states with twisted boundary conditions is easily computed with $\epsilon$-expansion due to the absence of zero modes, and the calculations done in the Appendices have included arbitrarily twisted boundary conditions to allow application to O($N)^{\ast}$ models.

Acknowledgements

S.W. and S.S. thank S. Chatterjee and A. Thomson for numerous helpful conversations. M.S., L.-P.H. and A.M.L acknowledge support through the Austrian Science Fund SFB FoQus (F-4018). The computational results presented have been achieved in part using the Vienna Scientific Cluster (VSC) and the HPC infrastructure LEO of the University of Innsbruck. S.W. and S.S. were supported by the NSF under Grant DMR-1360789, and MURI grant W911NF-14-1-0003 from ARO. Research at Perimeter Institute is supported by the Government of Canada through Industry Canada and by the Province of Ontario through the Ministry of Research and Innovation. SS also acknowledges support from Cenovus Energy at Perimeter Institute.

Appendix A Infinite volume computations

Here we recall some properties of the Wilson-Fisher fixed point in an infinite volume. We will need the resulting expressions to renormalize the theory in a finite volume, and we relate the couplings to universal observables of the infinite-volume theory. The computations are standard Zinn-Justin (2002), so we will be brief. Here we express the theory as a Euclidean action,

[TABLE]

where $d=3-\epsilon$, but now we integrate over infinite volume. At one-loop, we need to renormalize the couplings $s$, $u$, and the energy density $\Lambda$ (there is no field renormalization until two-loop). To this end, we introduce the renormalized couplings

[TABLE]

where $S_{d}=2/(\Gamma(d/2)(4\pi)^{d/2})$ is a convenient factor, and $\mu$ is an arbitrary renormalization scale. The renormalized coupling $s$ has been defined so that $s=0$ at criticality by definition. We renormalize the theory using a modified minimal subtraction scheme as detailed in Ref. Zinn-Justin (2002), where the renormalization constants are given by

[TABLE]

and the critical coupling is $s_{c}=0$. The Wilson-Fisher fixed point is obtained when the couplings take the values

[TABLE]

We are also interested in the vacuum energy in the vicinity of the Wilson-Fisher fixed point. To this end, we have included the additive constant $\Lambda$ which we choose to make the ground state energy density finite in an infinite volume. This constant also depends on the renormalization scale, and is given by

[TABLE]

Our choice of regularization leads to the vanishing of the ground state energy density at $L=\infty$ and $s=s_{c}$, where the system has full conformal invariance.

We note that computations will involve the arbitrary energy scale $\mu$. This dependence can always be eliminated in favor of physical quantities. As an example, we can consider the exact energy gap in an infinite volume when $s>0$. The inverse propagator can be written as a function of the Euclidean momentum $p=(\omega,k)$ as

[TABLE]

Then the exact energy gap $m$ is given by

[TABLE]

At one loop this gives

[TABLE]

where we have taken $g=g^{\ast}$. We can then always use this to rewrite dependence on $s$ and $\mu$ in terms of the physical parameter $m$ by inverting this expression,

[TABLE]

In an infinite volume, we can calculate the vacuum energy density. Setting $g=g^{\ast}$, this is given by

[TABLE]

where $E_{\gamma}$ is the Euler-Mascheroni constant and $s>0$.

We note that our finite-volume calculations must be analytic through the critical point $s=0$, and the system remains disordered for $s<0$. Therefore, any singularities or branch cuts present in these expressions must cancel out in final results. Therefore we prefer to give our expressions in terms of $s$ and $\mu$ rather than $m$ and $E_{0}$, which are both non-analytic for $s=0$, and undefined for $s<0$.

Appendix B Derivation of the Bloch effective Hamiltonian

Here, we give a derivation of Eq. (14), the effective Hamiltonian from Bloch’s perturbation theory Bloch (1958). For the relation between this approach and other effective Hamiltonian approaches to perturbation theory, we direct the reader to Ref. Klein (1974).

We begin by considering a degenerate subspace of the unperturbed spectrum, $\Omega_{0}=\mathrm{Span}\{|\epsilon_{0}\rangle\}$, where the states $|\epsilon_{0}\rangle$ satisfy

[TABLE]

These states are split into distinct energies by the exact Hamiltonian $H=H_{0}+V$,

[TABLE]

Here, we have defined the orthonormal basis of states $|\alpha\rangle$ which diagonalize the exact Hamiltonian and reduce to the degenerate manifold $|\epsilon_{0}\rangle$ in the $V=0$ limit. We will define the space spanned by these states by $\Omega=\mathrm{Span}\{|\alpha\rangle\}$. If the perturbation is small, there should be a one-to-one correspondence between the spaces $\Omega_{0}$ and $\Omega$. We will assume the latter fact in what follows.

Let $P_{0}$ be the projection operator onto the space $\Omega_{0}$, and define the states

[TABLE]

The set of states $|\alpha_{0}\rangle$ are a particular basis spanning $\Omega_{0}$, although this basis is not orthonormal in general. We now claim that we can define a linear operator $U$ such that

[TABLE]

Once the operator $U$ is found, we can construct an effective Hamiltonian which acts on the unperturbed subspace but gives the exact energy spectrum,

[TABLE]

This Hamiltonian acts on the unperturbed subspace, but gives the exact spectrum of the interacting Hamiltonian. We can also obtain the exact eigenstates by $|\alpha\rangle=U|\alpha_{0}\rangle$.

We now make a few clarifying comments on the above operators and states. We note that $\Omega_{0}$ and $\Omega$ represent subspaces within the same Hilbert space, and that these subspaces overlap by assumption. Since the linearly independent basis $|\alpha_{0}\rangle$ defined via (84) is not necessarily orthogonal, the operator $U$ is not unitary, so the standard intuition on changes of basis does not apply here. For example, since we can decompose any state $|\phi\rangle$ in the Hilbert space as

[TABLE]

it follows that

[TABLE]

so

[TABLE]

is an exact operator identity. In particular, this implies $U|\alpha\rangle=|\alpha\rangle$. Furthermore, the effective Hamiltonian defined in Eq. (86) is not necessarily Hermitian. Thus, even though the vectors $|\alpha_{0}\rangle$ are eigenvectors of $H_{\mathrm{eff}}$, it does not follow that the effective Hamiltonian is diagonal in the $|\alpha_{0}\rangle$ basis. This can also be inferred from the fact that the $|\alpha_{0}\rangle$ may not be orthogonal. Below, we will define a similarity transform which does allow the definition of a Hermitian effective Hamiltonian from Eq. (86).

We now find an explicit expression for $U$ for a completely general Hamiltonian. Starting from Schrödinger’s equation,

[TABLE]

Since $UP_{0}=U$,

[TABLE]

for all $|\alpha\rangle\in\Omega$. Furthermore, for any state $|\alpha\rangle\in\Omega$, we have $U|\alpha\rangle=|\alpha\rangle$ and $UH_{0}|\alpha\rangle=\epsilon_{0}|\alpha\rangle$, from which we have the general operator relation on the entire Hilbert space:

[TABLE]

We can invert the left-hand side of this equation by using the fact that the operator $\left(\epsilon_{0}-H_{0}\right)$ satisfies $\left(\epsilon_{0}-H_{0}\right)P_{0}=0$ and has a well-defined inverse in the space orthogonal to $\Omega_{0}$. This gives the implicit equation

[TABLE]

This allows an expansion in powers of $V$. Up to third order in $V$, this is given by

[TABLE]

For higher order expressions, see Refs. Bloch (1958); Lüscher (1983). As seen in Section III.2, the leading non-trivial two-loop diagram is of order $\epsilon^{5/3}$. We can now check explicitly by putting $V$ from Eq. (10) that the last term in Eq. (94) is of order $\epsilon^{5/3}$ or higher, justifying our truncation of the effective Hamiltonian in the main body of our paper.

We note that the effective Hamiltonian defined above will not be Hermitian in general without a redefinition, although we will not encounter this problem at one-loop. In general, we have

[TABLE]

where $P$ is the projection operator onto the space $\Omega$. Then the redefined Hermitian Hamiltonian

[TABLE]

acts on $\Omega_{0}$ and has the same spectrum as $H_{\mathrm{eff}}$.

Appendix C Loops sums

In this appendix we give the calculation of the relevant loop diagrams on a torus in fractional dimensions. The torus is parametrized by complex coordinates, $x=x_{1}+ix_{2}$, with two complex periods $\omega_{1}$ and $\omega_{2}$, see Fig. 1. We define the modular parameter, $\tau=\omega_{2}/\omega_{1}=\tau_{1}+i\tau_{2}$, and the length scale $L\equiv|\omega_{1}|$. The area of the torus is given by $\mathcal{A}=\tau_{2}L^{2}$.

In this geometry, the basis vectors of the dual lattice are

[TABLE]

We consider the eigenvalues of the Laplacian on the torus. In the interest of generality, we will allow twisted boundary conditions along the two cycles of the torus.

[TABLE]

The parameters $a_{1}$ and $a_{2}$ take values in the range $[0,1)$. For $a_{1,2}$ not equal to either an $1$ or $1/2$, the fields $\phi_{\alpha}$ are actually complex and our symmetry breaks down from O($N$) to $SU(N)\times U(1)$. In the main text of this paper we always take $a_{1}=a_{2}=0$, but the calculation of loop sums on the torus with twisted boundary conditions has found applications in related work Schuler et al. (2016); Whitsitt and Sachdev (2016); Whitsitt et al. (2016), so we give general results here.

The eigenvalues of the Laplacian are

[TABLE]

and a general one-loop diagram will take the form

[TABLE]

where $\gamma=\tau_{2}^{2}L^{2}s/4\pi^{2}$, and the primed summation indicates that the $n=m=0$ term is omitted in the fully periodic case $(a_{1},a_{2})=(0,0)$.

We now generalize this sum to arbitrary dimension. This is done by promoting the two-vector $(n+a_{1},m+a_{2})$ to a $d$-dimensional vector where the first $d/2$ components are $n+a_{1}$ and the last $d/2$ components are $m+a_{2}$. Then in (100) the sums are taken over $n,m\in\mathbb{Z}^{d/2}$. We write the sums as

[TABLE]

The summand is rewritten using the identity

[TABLE]

giving

[TABLE]

We can now write the sum in terms of the two-dimensional Riemann theta function, defined as

[TABLE]

where $\Omega$ is a $2\times 2$ matrix and $\mathbf{u}$ is a two-dimensional vector. Then

[TABLE]

where

[TABLE]

The function $S_{\nu}^{(d)}(s,\tau)$ has a divergence for small $\lambda$ whenever $\nu<d/2$. We evaluate the sums by dimensional regularization: we separate out the divergent parts and and evaluate them in the convergent regime $\nu>d/2$, and then analytically continue them to the physical values of $\nu$ and $d$ of interest, taking care to renormalize any poles in $\epsilon$ which arise. For our purposes, it is also crucial that we obtain final expressions which are regular at $s=0$ and remain finite for $s<0$, since the finite-volume theory should be analytic through the critical point.

We proceed by splitting the integral in Eq (105) into two parts, $\int_{0}^{\infty}=\int_{0}^{1}+\int_{1}^{\infty}$, and studying the divergent part. Using the mathematical identity

[TABLE]

we write the lower portion of the integral as

[TABLE]

where we define the vector $\mathbf{v}_{2}=-i(a_{1},a_{2})$. For $d=3$ and $\nu\geq-1/2$, the UV divergences are entirely contained in the last line of Eq (108). The very last term can be integrated in its convergent regime, obtaining

[TABLE]

This expression was obtained by evaluating the integral for $\nu>0$ and $\gamma>0$, obtaining an expression involving incomplete gamma functions, and then expressing these functions as a power series. The final expression is a single-valued analytic function for all $\gamma$ and $\nu$ with no singularities Abramowitz and Stegun (1964), and evaluating this series numerically is trivial.

When $\nu>d/2$ we can evaluate

[TABLE]

Since we are expanding around $d=3$, the three terms will contribute poles for $\nu=3/2$, $\nu=1/2$, and $\nu=-1/2$ respectively. We will see that these three poles are related to the renormalization of $u$, $s_{0}$, and $\mathcal{E}_{0}$.

To summarize our results so far, we have written the loop sum as

[TABLE]

[TABLE]

where $\gamma=\tau_{2}^{2}L^{2}s/4\pi^{2}$, the Riemann theta function $\Theta$ is defined in (104), and

[TABLE]

From Eq. (112) we define finite functions. They are given in terms of the $S_{\nu}^{(d)}$, but with their poles around $d=3$ removed.

[TABLE]

The three functions $f_{\nu}^{(d)}$ have a regular power series about $s=0$ and $\epsilon=0$, and the extra factors are defined to simplify expressions in our renormalization scheme. We note that these functions are all dependent on $\mu$, but the first two functions are independent of $\mu$ at $s=0$, while the third function will be exchanged for a $\mu$-independent function below. We note the identities

[TABLE]

The first identity will be useful for renormalizing the bare mass, while the other two are a consequence of modular invariance.

Appendix D Renormalization of the effective Hamiltonian

Here we consider the renormalization of physical quantities. To this end, we calculate the effective Hamiltonian for the splitting of the Fock vacuum, $h_{k=0}$, which contains most of the divergences which need to be considered.

We calculate each individual expression in Eq. (30). The constant term is

[TABLE]

By expanding out all of the terms, including the $g$-dependent $Z_{2}$ factor using Eq. (74), and using the definition of $\Lambda$ from Eq. (76), we find that all poles in $\epsilon$ cancel, as well as all factors of $S_{4-\epsilon}$. Once the poles cancel, we set $g=g^{\ast}$, obtaining

[TABLE]

The coefficient of $\varphi_{\alpha}^{2}$ is

[TABLE]

The poles in $\epsilon$ and factors of $S_{4-\epsilon}$ cancel, and after setting $g=g^{\ast}$ we obtain

[TABLE]

Finally, the divergence in the quartic term cancels similarly, and after setting $g=g^{\ast}$ we find

[TABLE]

We note all of these expressions appear to depend on the arbitrary scale $\mu$. This dependence actually drops out of $h^{(4)}_{k=0}$ to this order in $\epsilon$, and it drops out of all quantities at the critical point $s=0$. We can write

[TABLE]

where $f^{(3)}_{3/2}(\tau,s)$ is $\mu$-independent. The fact that all $\mu$ dependence vanishes at $g=g^{\ast}$ and $s=0$ is a manifestation of the scale invariance of the critical theory.

Combining the above results, the effective Hamiltonian takes the form given in Eq. (31):

[TABLE]

To the order required, the special functions are

[TABLE]

[TABLE]

[TABLE]

The function $f_{-1/2}^{(3-\epsilon)}$ should be expanded to first order in $\epsilon$. It is possible to exchange the parameters $s$ and $\mu$ for the infinite volume gap and ground state energies, but the latter are not analytic through the critical point so we keep the $\mu$ dependence in our final expressions.

Appendix E Correspondence with the large-$N$ expansion

In a recent paper, two of us have computed the spectrum of the O($N$) model at $N=\infty$ Whitsitt and Sachdev (2016). In this approach, one starts from a saddle-point of the theory which already has an energy gap on the torus, so the theory is IR safe and the zero mode plays no special role in the expansion. In this Appendix we take the small-$\epsilon$ limit of our large-$N$ results and compare them to the large-$N$ limit of the $\epsilon$-expansion. We find exact agreement where possible in both methods. Here we limit ourselves to $s=0$.

In the large-$N$ expansion, one begins by solving the gap equation on the torus, which can be written

[TABLE]

We can evaluate this using the methods of Appendix C, where dimensional regularization sets the integral on the right-hand side to zero, obtaining from Eq. (115)

[TABLE]

As $\epsilon\rightarrow 0$, the function $f^{(3)}_{1/2}$ is regular, while the last term in Eq. (127) diverges. This requires the first term to diverge in the same way, from which one already sees that to leading order $\Delta\sim\epsilon^{1/3}$. Continuing this process, one can explicitly solve Eq (127) perturbatively in $\epsilon$, making use of identities derived in Appendix C. The gap at $N=\infty$ is given by

[TABLE]

where the functions are identical to those defined in Eqs. (123-125), and the $\mu$-dependence has dropped out.

Once the gap equation has been solved, the Hamiltonian of the theory at $N=\infty$ is given by

[TABLE]

Here,

[TABLE]

which should be evaluated order-by-order in $\epsilon$ using the techniques in Appendix C. The operator $P$ projects onto all states which contain no O($N$) singlets. Finally, the energies of the singlet states $E_{n}(k)$ are the solutions to the equation

[TABLE]

where

[TABLE]

We now consider a subset of the above spectrum; specifically, the states created by

[TABLE]

This Hamiltonian creates two kinds of zero-momentum “particles,” with masses $\Delta$ and $E_{1}(0)$ respectively. The particle with mass $E_{1}(0)$ transforms as an O($N$) singlet, while the states with $\ell$ mass-$\Delta$ particles are in the $\ell$th symmetric traceless tensor representation of O($N$). We now argue that $H_{\mathrm{eff},k=0}$ is precisely the $N=\infty$ limit of the Fock vacuum Hamiltonian in the $\epsilon$-expansion.

We saw in Section IV that the solutions to the Fock vacuum Hamiltonian satisfy the equation

[TABLE]

Recall that this describes states in the $\ell$th symmetric traceless tensor representation of O($N$). We now use the large-$N$ expansion in quantum mechanics. The idea is that the centrifugal term acts as an effective mass at large-$N$, resulting in a harmonic well at the stationary point of the effective radial potential provided it is well-behaved. For a review of this expansion, see Ref. Chatterjee (1990), which gives an explicit formula (3.2.20) for the spectrum to the first few orders in $1/N$. At $N=\infty$ we find the spectrum

[TABLE]

where $\mathcal{E}_{0}$ and $\Delta$ agree exactly with their expressions calculated in the large-$N$ expression Eq. (133) to order $\epsilon^{4/3}$. Furthermore, while we cannot compare $E_{0}$ directly with the first zero of $\Pi(0,\omega)$, we have evaluated it numerically for small values of $\epsilon$ and found very good agreement. Finally, the irreducible representations of the states under the O($N$) symmetry agrees exactly.

We mention that the $1/N$ expansion in quantum mechanics is an easy way to obtain $1/N$ corrections to the Fock vacuum Hamiltonian compared to the field-theoretic methods in Ref. Whitsitt and Sachdev (2016). One may also attempt to calculate the spectrum for the $k>0$ effective Hamiltonians using the large-$N$ expansion, which has been successfully applied to non-isotropic Hamiltonians in atomic and molecular physics Chatterjee (1990).

Appendix F Strong-coupling expansion of isotropic quartic oscillators

In this appendix we give details of the numerical calculation of the spectrum of the isotropic quartic oscillator

[TABLE]

in the strong-coupling limit, finding the coefficients

[TABLE]

We tabulate the values of $c_{n,\ell,m}$ which we have calculated in Tables 4-6.

We begin by solving Eq. (136) numerically for $r=0$. We first fix the asymptotic behavior by writing

[TABLE]

The exponential factor takes into account the large-$\rho$ behavior implied by Eq. (136). We then use a shooting method, using known boundary conditions on the wave function at $\rho=0$ (including an arbitrary normalization) to compute $\psi_{n,\ell}$ in Mathematica using DSolve for variable values of the energy until we identify an eigenstate.

Once we find the energy to sufficient accuracy, the function $\psi_{n,\ell}$ will not change much for smaller $\rho$, but will always blow up after some value of $\rho$. However, the actual wave function $R_{n,\ell}$ is exponentially suppressed, so we only need to obtain $\psi_{n,\ell}$ accurately for small values of $\rho$ to obtain an accurate wave function.

Once we obtain a numerically accurate energy and wave function for $r=0$, we use logarithmic perturbation theory to compute the expansion in $r$. This has the benefit of only needing the unperturbed energy and wave function, whereas the standard Rayleigh-Schrödinger expansion requires knowledge of many excited states to get accurate values for the coefficients. The starting point for logarithmic perturbation theory is to write the wave function as

[TABLE]

where $\rho_{i}$ are the nodes of $R_{n,\ell}$, and then rewrite the eigenvalue equation as an equation for $G$. Then we assume an expansion in $r$:

[TABLE]

Inserting these definitions into our eigenvalue equation, the resulting differential equation is linear order-by-order in perturbation theory, so the energies are given in closed form in terms of integrals only involving the unperturbed functions $G_{0}(\rho)$, $E_{n,\ell,0}$, and $\rho_{i,0}$. For explicit details, we refer the reader to References Au and Aharonov (1979); Aharonov and Au (1979). We give the coefficients we have calculated in Tables 4-6.

Appendix G Numerical Simulations

To compute the critical torus energy spectrum numerically with exact diagonalization (ED) we consider explicit lattice models with a critical point known to be in the respective universality class. Typical Hamiltonians for such models with critical points in the O(2) and O(3) universality classes were given in Eqs. (68) - (71).

In a first step, we calculate the low-energy spectrum of a given Hamiltonian on finite periodic clusters with up to $N=36$ lattice sites at the respective quantum critical point. The spectrum can be divided into $S^{z}$ symmetry sectors combined with the irreducible representations of the lattice space-group symmetry. The critical low-energy spectrum collapses as $1/\sqrt{N}=1/L$ with the linear system size $L$. We multiply the spectrum by $L$ to get rid of this scaling. Finally, we extrapolate the scaled spectrum for the different system sizes linearly in $1/N$ to the thermodynamic limit to obtain the universal critical torus energy spectrum. This extrapolation approach was successfully used and corroborated by large-scale quantum Monte Carlo simulations in Ref. Schuler et al. (2016) for the Ising CFT.

In Fig. 10 we show the torus spectrum for the considered finite square clusters with $N$ lattice sites for the model Eq. (68). The different panels show results for different values of $J_{z}$. After the linear extrapolation in $1/N$ to the thermodynamic limit the spectra for the different $J_{z}$ agree to high precision with each other, demonstrating its universal nature. As we cannot easily determine the speed of light for the models with ED, we normalize all levels such that the extrapolated energy gap to the first excited state with $\ell=0$ is set to one, $\Delta_{\ell=0}\equiv\Delta_{\varepsilon_{T}}=1$. After this, we have finally obtained the critical torus energy spectrum from numerics, which is displayed in Fig. 7.

The same procedure is used to obtain the critical torus spectra for the other Hamiltonians (not shown). There are, however, some subtleties to consider. First, for the O(3) models $\ell=S_{\text{tot}}$ is usually not resolved in ED. Nevertheless, $S_{\text{tot}}$ can usually be obtained by detecting levels with the same energy with quantum numbers $-S_{\text{tot}}\leq S^{z}_{\text{tot}}\leq S_{\text{tot}}$. Furthermore, for non-trivial magnetically ordered phases on one side of the critical point, like the orderered state in model Eq. (70), all space-group symmetry sectors in the corresponding tower of states appear as levels in the $\kappa=0$ sector. For the model Eq. (70) on the square lattice this means that the levels in the $\kappa=0$ sector consist of levels with momenta $\mathbf{k}=(0,0,0)$ and the Bragg peak momentum $\mathbf{k}=(0,0,\pi)$. Subsequently, the $\kappa=1$ levels are composed of the momenta closest to $(0,0,0)$ as well as the momenta closest to $(0,0,\pi)$ etc. In QMC simulations the antiferromagnetic version of this model with $J>0$ is usually considered. Then the $\kappa=0$ levels would include levels with $\mathbf{k}=(0,0,0)$ and $\mathbf{k}=(\pi,\pi,\pi)$.

In the following tables 7 - 10 we list the numerical values for the critical torus spectra for the O(2) and O(3) CFTs on both, square and triangular lattice geometries. Numerical values for the Ising CFT can be found in the Supplemental Material of Ref. Schuler et al. (2016).