Conformal dimensions via large charge expansion

We construct an efficient Monte Carlo algorithm that overcomes the severe signal-to-noise ratio problems and helps us to accurately compute the conformal dimensions of large-$Q$ fields at the Wilson-Fisher fixed point in the $O(2)$ universality class. Using it we verify a recent proposal that conformal dimensions of strongly coupled conformal field theories with a global $U(1)$ charge can be obtained via a series expansion in the inverse charge $1/Q$. We find that the conformal dimensions of the lowest operator with a fixed charge $Q$ are almost entirely determined by the first few terms in the series.

Conformal field theories (CFTs) occupy a central place in our understanding of modern physics. They describe critical phenomena in condensed matter physics and statistical models [1,2], quantum gravity via the AdS/CFT correspondence [3] and can be found at fixed points of renormalization group flows [4][5][6][7]. They are uniquely described by a set of dimensionless numbers (the CFT data), i.e. conformal dimensions and OPE coefficients associated with the primary fields of the theory. Since they are typically strongly coupled and lack a characteristic scale, it is often difficult to compute the conformal dimensions analytically. Still, a number of sophisticated techniques have been developed to deal with this challenge, both perturbatively (e.g. 4 − expansion, fixed-dimension expansion, large-N ; see [8] for a review) and non-perturbatively (e.g. bootstrap [9]). In some cases Monte Carlo techniques offer a reliable numerical approach for computing the conformal dimensions [10,11].
Energies of low-lying states also capture universal features of a 2 + 1 dimensional CFT when the theory is studied on a space-time manifold R × Σ (see [12,13] for some recent work in this direction). For example conformal dimensions D of operators on R 3 are related to the energies E S 2 of states living on a two-sphere of radius r 0 through the relation D = r 0 E S 2 [14,15]. This relation, known as the state-operator correspondence, is a consequence of the fact that R 3 is conformally equivalent to R × S 2 (r 0 ). Recently, such a connection has been used in CFTs with global U (1) charges to show that the conformal dimension D(Q) of the lowest operator with fixed U (1) charge Q can be expanded in inverse powers of the charge density on a unit sphere Q/4π [16,17] (see also [18][19][20][21][22] for related work): where c 0 ≈ −0.094 [23] and the other coefficients only depend on the universality class. While a simple dimensional analysis allows one to predict the leading large-Q behavior, it is a priori unclear if a power series can capture the sub-leading corrections. The recent work argues that by separating the theory into sectors of fixed charge Q one can construct an effective field theory (EFT) in each sector, which can be used to compute the energies as a power series in 1/Q. Through the state-operator correspondence one can then obtain the series expansion of the conformal dimension D(Q) and relate directly the coefficients in Eq. (1) to the coefficients in the energy expansion. In this work we make significant progress in establishing that Eq. (1) is an excellent description of the WF fixed point in the O(2) universality class. In order to achieve this we overcome severe signal-to-noise ratio problems in Monte Carlo methods that usually hinder calculations of D(Q) for large values of Q. Our new approach allows us to determine the corresponding universal coefficients arXiv:1707.00711v3 [hep-lat] 14 Jan 2018 c 3 2 = 1.195(10), c 1 2 = 0.075(10) for the first time. In Fig. 1 we demonstrate that the measured values of D(Q) using our Monte Carlo approach are excellently described by Eq. (1). Surprisingly, we find that the large-Q expansion with the first three coefficients matches the conformal dimensions even at Q = 1 within a few percent. Thus, our work demonstrates that at least in some class of models, the large Q expansion is similar to the epsilon expansion in the fact that the strongly-coupled conformal fixed point can be described within a simple perturbative framework. The coefficients of the expansion could still be difficult to compute analytically, but perhaps bootstrap techniques could be developed for it [24].
To understand the origin of Eq. (1), consider the conformal field theory describing the Wilson-Fisher fixed point in the three-dimensional O(2) universality class. In a fixed charge sector Q, the charge density introduces the mass scale Q/V in the theory and hence for momentum scales p Q/V the physics is described by a Goldstone field χ that is controlled by an approximately scale-invariant Lagrangian [16,18] (see also [25] for a related approach to effective descriptions of non-relativistic CFTs): where R is the scalar curvature of the manifold R × Σ. Thus, we learn that in the large-Q limit only the two parameters k 3 2 and k 1 2 that appear in Eq. (2) play an important role and all other terms are suppressed [17,19]. Since the charge is non-zero, the action is only meaningful away from χ = 0 and is to be expanded around the fixedcharge homogeneous classical solution χ = µt. Using the effective quantum Hamiltonian arising from the effective Lagrangian Eq. (2), one can show that the total energy of the system is given by where the first two terms are related to the couplings in the effective Lagrangian Eq. (2) through the relations The higher order terms in the expansion are related to higher dimensional operators in Eq. ((2)) and quantum corrections. The last term q Σ arises due to quantum fluctuations that can be computed exactly for simple manifolds. For the sphere (R = 2/r 2 0 ) one finds q S 2 = c 0 /r 0 where c 0 ≈ −0.094 [23], while for the torus (R = 0) it is q T 2 = c 0 /L with c 0 ≈ −0.508 [26].
By choosing Σ = S 2 and using the state operator correspondence one can now easily derive Eq. (1). It is interesting to note that the coefficients c 3 2 , c 1 2 in Eqs. (1),(3) are related to the low-energy constants k 3 2 and k 1 2 of the effective Lagrangian in Eq. (2). Indeed, these low-energy constants are independent of the manifold chosen and depend only on the CFT. Assuming the manifold is the torus we predict that Note that every term in the energy expansion is a dimensionless function of three variables: a coefficient in the D(Q) expansion, a geometrical term from the manifold, and a power of V /Q. The motivation of our current work is to compute D(Q) and E Σ (Q) in the classical O(2) sigma model on a torus and verify the expansions in Eq. (1) and (3) and the connections between them. We accomplish this by regularizing the classical O(2) sigma model on a cubic lattice with lattice spacing a and use Monte Carlo methods to perform the calculations. The model is defined by phases, exp(iθ x ) on each three-dimensional lattice site x = (x 1 a, x 2 a, x 3 a) and the nearest neighbor action Hereαa denotes the three unit lattice vectors, and β is the coupling of the model. The physics of the Wilson-Fisher fixed point can be studied by tuning the coupling to its critical value (β c = 0.4541652 [27][28][29]), where a second-order phase transition separates the symmetric phase (β < β c ) and the spontaneously broken superfluid phase (β > β c ). Universality implies that details our specific model should be irrelevant in the limit a → 0 which is naturally reached by studying large lattices at β c .
Configurations that contribute to the partition function of the lattice model at the critical point can be efficiently generated by both the Wolff cluster algorithm [30] and the worm algorithm based on the worldline representation [31]. However, in order to compute the conformal dimension D(Q) in R 3 we need to compute the two-point correlation function C Q (r) of charge Q fields on a large lattice of size L, which is expected to decay as a power law for large separations r L at the critical point: Fitting the data to this form we can in principle extract D(Q) and thus verify Eq. (1). Note that for Q = 1, it reduces to the standard 2-point correlation function, which is used to extract the critical exponent η through the relation C 1 (r) = G(r) ∝ 1/r d−2+η . For Q = 2, 3, 4, the corresponding conformal dimensions have also been computed earlier using different methods, and the results are summarized in Table I [33] for Q = 3 and in [34] for Q = 4), previous MC results are in column four [35], and bootstrap results are in column five [36]. With the Wolff cluster algorithm it is difficult to average numbers of order one to compute a small value of C Q (r) at large separations. In contrast, in the worm algorithm, it is difficult to correctly build the worldline configurations that contribute to the correlation function in the presence of charged sources Q and −Q separated by a large distance. In this case the severe signal-to-noise ratio problem emerges as an overlap problem between the vacuum ensemble and the one containing the sources.
We can now verify if the conformal dimensions in Table II are consistent with Eq. (1). For this purpose we perform a combined fit of our data for the difference ∆(Q) and D(Q) assuming that c 3 2 , c 1 2 , c − 1 2 are non-zero and c 0 = −0.094 as expected. Taking into account various systematic errors from fitting procedures we estimate c 3 2 = 1.195(10), c 1 2 = 0.075(10) and c − 1 2 = 0.0002(5). The raw data for ∆(Q) are shown in Fig. 3, and further technical details are discussed in the Supplementary Material. We also show a comparison with the prediction obtained by just keeping the first three leading terms of the expansion in Eq. (1). As the figure shows, this prediction works even at small values of Q but is off only by a few percent at Q = 1.
Next we explore if we can connect our above calculations of c 3 2 and c 1 2 with the ones appearing in Eq. ((3)) for the expansion of the energy on a torus. Lattice calculations naturally lead to a torus geometry in the continuum limit if we keep the physical length L fixed while taking the number of lattice points in each direction, L/a, to infinity. The lattice spacing a itself is defined by setting the lattice energy E L (Q) to be equal to the continuum energy E T 2 (Q) on the torus of side L as the continuum limit is taken. On the lattice we measure the energy in terms of the dimensionless number E L (Q)a t as a function of L/a, where a t is the temporal lattice spacing. Although the lattice calculation at a fixed L/a and charge Q breaks the symmetry between space and time, in the continuum limit (i.e., L/a → ∞), we expect a t /a → 1 due to the cubic symmetry of the lattice action, Eq. (5). Thus on the lattice with a fixed L/a we can replace Eq. (3) by the formula, such that in the continuum limit (i.e., large lattices) we expect this equation to turn into Eq. (3) on the torus, Unfortunately, computing E L (Q) on the lattice is subtle due to the additive renormalization of lattice energies [37][38][39]. Hence in this work we use the techniques discussed in [31] to compute energy differences ∆E L (Q) = (E L (Q) − E L (Q − 1)) at a fixed lattice size L/a. The idea is to couple a chemical potential µa t to the conserved charge Q and extract the energy differences between ground states with charges Q and Q − 1 by tuning the chemical potential. At the critical chemical potential µ (Q−1) c = ∆E L (Q) the average charge of the system jumps from Q − 1 to Q. Further details on the method can be found in the supplementary material.
Unfortunately, our method is efficient only on small lattice sizes L/a, limiting the range of Q's that can be used. Remember that the EFT description is valid only when 1 (L/a)/ √ Q L/a. Further, small lattices also imply larger lattice spacing errors which means c L   Table III as a function of Q along with a plot of the fitted function (solid line) with the fit parameters c L Our data for L/a = 8 and 10 fits well to Eq. (8) with a χ 2 /DOF ≈ 1, as long as we restrict the fits to the range 5 ≤ Q ≤ 13 and 10 ≤ Q ≤ 18 respectively. The fit result for L/a = 10 is shown in Fig. 4. The fits give c L is only about 3% off from 1.195(10) extracted earlier using conformal dimensions. This error seems reasonable given our small lattices (see the supplementary material for further discussion of this point). The coefficient c L 1 2 is small, perhaps related to the fact that it is expected to be zero in the continuum. The value of c − 1 2 is consistent with being non-zero, although it was recently argued that the classical contribution to the energy at this order is expected to vanish in the continuum limit [19]. Contributions from quantum corrections were not computed and could in principle be non-zero. Finally, as discussed in the supplementary material, all of our results are consistent with those obtained recently [13]. As explained in [31], we can construct a worm algorithm to compute quantities in statistical mechanics, starting with the classical action We first use the identity on each bond (x, α) to integrate out the field variables (θ x ) from the partition function and express it in terms of a configuration of integer bond variables k x,α , each of which is an integer valued worldline variable denoting the charge flowing on the bond (x, α), from x to (x +α). The function I k (β) is the modified Bessel function of the first kind. In terms of worldline configurations [k], the partition function looks like (12) We can use the standard worm algorithm to update worldline configurations [k] as follows: 1. We pick a random site on the lattice, x = x h and refer to it as the head site. At the beginning we also define the tail site to be the same as the head site, x t = x h .
3. Let k be the current on the bond joining x t and x t +α. If α > 0, we update the forward current from k to k + 1 with probability I k+1 (β)/I k (β). If this update is accepted, then we change the tail site from x t → x t +α. If α < 0, then with probability I k−1 (β)/I k (β), we update k to k − 1 and set x t → x t +α. When the updates are not accepted, x t remains the same.
4. If the tail site x t reaches the head site x h the worm update ends. Otherwise, we go back to step 2.
In the above update, when we start we introduce a single positive charge at x h and a single negative charge at x t . During the worm update the tail site x t moves around the lattice until it comes back and annihilates with the We can easily adapt the worm algorithm to measure the two point correlation function C Q (r) = e iQθ0 e −iQθr . We simply introduce charge Q at the head site and a charge −Q at the tail site by updating the bonds by changing k to k ± Q as we move. Whenever the tail reaches the site x t = x h + r, we get a contribution to C Q (r). Combining this with the standard worm algorithm we can in principle get an ergodic algorithm. However, such an algorithm leads to a severe signal-to-noise problem, since the updates of k to k ± Q is extremely inefficient especially for large Q.
To solve this problem we focus on the ratio of correlation functions where we define with Q x = Q(δ x,o − δ x,r ). In other words Z Q,r is a partition function with a charge Q fixed at x = 0 and −Q at x = r. Note that we can also write the above ratio slightly differently as where the expectation value on the right is now taken in the distribution of [k] according to the partition function Z Q,r . In this case the standard worm algorithm can be used to update the configurations [k] associated with Z Q,r , which we refer to as the target ensemble. In Fig. 5 (left) we show an illustration of a configuration in the target ensemble with a worm loop that updates the background configuration. In addition to the standard algorithm, we construct a special measurement algorithm where the worm update starts and ends at x h = 0. If during this measurement update the tail site touches reaches the site x t = r we count it as a contribution to the ratio R(r). More concretely, 1. We start with a configuration which has Q charges at site (0, 0, 0) and −Q charge at site (r, 0, 0), and initialize a counter c = 0.
2. We begin a worm update with the head at x h = (0, 0, 0). We also define the tail site as x t = x h .
3. We pick one of the 2d neighbors of x t , and propagate the tail site using the standard worm update described earlier. 4. Whenever x t = (r, 0, 0), we increment the counter c = c + 1. This implies that the configuration generated contributes to the ratio R(r).

5.
When the tail site returns to the head site the update ends.
An illustration of this measurement update is shown in Fig. 5 (right), where the worm starts from x h = 0 and passes through x t = r before closing. R(r) is computed as an average of the value of c measured for each measurement worm algorithm. The actual algorithm involves several standard worm updates and an equal number of measurement updates. For our simulations, we work at several finite cubic volumes with linear size L and compute the correlation ratios with r = L/2. In this calculation, the flux from charge Q can reach the sink either directly through the bulk or through the boundary. We can define a winding number w for each configuration [k] as the number of flux lines that reach the sink through the boundary. Then the remaining Q − w flux lines reach through the bulk. We have discovered that the value of R(r) is sensitive to w. In our work we sample all possible windings distributed according to Z Q,r . On an average we see that w = Q/2. The standard worm algorithm is able to change the winding number w quite efficiently.

DETAILS OF FITTING AND ERROR ANALYSIS
In the main text, we have presented evidence for the validity of the large charge analysis developed in [16,17]. In this Supplementary we first explain the analysis details to compare the theory with the Monte Carlo data, as well as model (in)-dependence of our results. We also compare our results with some related results in the literature [13] for the U(1) global charge in (2 + 1)−dimensions.

Conformal Dimensions
The main text, and the first supplementary section describes rather extensively the new method developed to estimate the difference in the conformal dimensions, D(Q) − D(Q − 1) directly from the Monte Carlo simulations. Absolute conformal dimension D(Q) are obtained by using (with D(0) = 0): Since the individual terms, D(Q) − D(Q − 1), are obtained by fits to completely different sets of independent simulations, there is no correlation between the different sets and the errors can be added in quadrature to obtain the error on D(Q): A bootstrap analysis was also done to check the error-bars on these values. Fit methodology and model (in) dependence: As emphasized in the main text, we have tried to verify the validity of the effective theory developed in [16,17], which predicts the conformal dimensions D(Q) as in Eq. (1). However, from the Monte Carlo results we can directly measure only the difference D(Q) − D(Q − 1). Together with the input of D(Q = 0) = 0, we can reconstruct the values for D(Q). To verify the effective theory (EFT), we fixed the value of c 0 from the EFT and then fitted both the difference of the conformal dimensions, and the total conformal dimensions to their respective functional forms (Eq. (8)) to extract the low energy coefficients c 3 2 and c 1 2 simultaneously. The results of this fit are described in the main text. We note that the systematic power law expansion of the conformal dimensions in 1/Q is rather robust, it is possible to relax the input from the EFT for the value of the additive constant c 0 , and determine it completely from the numerical results by treating it as a fit parameter. For example, keeping the constants c 3 2 and c 1 2 fixed to their values as described in the methodology before, we have tried to fit c 0 from the D(Q) results. In this case, we obtained for c 0 the value -0.102(4), where the errorbar takes into account the different fit ranges, about a 2-sigma difference from the value predicted by the EFT. We also tried the different strategy of fitting all the coefficients c 3 2 , c 1 2 and c 0 from the data on D(Q). This yields: c 3 2 = 1.195(5), c 1 2 = 0.07(1) and c 0 = −0.09 (2). As one notes, these values are completely consistent with the results quoted in the main text, as well the input from the EFT. We note that these different strategies show that the assumed functional form for the conformal dimensions is rather robust, and gives confidence in the model independence of the results. Another important point in the analysis is that from the ansatz of Eq. (1), there is a contribution from the term c − 1 2 (4π/Q) 2 , which is however very suppressed, even with the accuracy of our data. This contaminates the extraction of either the term c 0 , or the estimate of c − 1 2 . In testing the different strategies, we have favored ranges with larger starting values, such that the latter term can be fixed to zero in the fit.

Energy
In this subsection, we explain the methodology for the energy computation, used in [31]. This is completely different method used for calculating the conformal dimensions and provide an independent way of checking the EFT. The basic idea is to couple a chemical potential µ to the conserved global U (1) charge Q. As the chemical potential is increased, the ground state is the one which has a total charge Q, since the chemical potential couples as H − µQ. The action for the XY-model in the presence of the chemical potential is Thanks to the worm algorithm, this model can be simulated with a chemical potential as demonstrated in [31]. It was verified that with increase in chemical potential there are level crossing phenomena between a state E (Q−1) with a charge Q − 1 and a state E Q . With increasing charge, the latter becomes the ground state of the system. Close to the critical chemical potential where the level crossing occurs, one can approximate the partition function as where E Q 0 is the ground state in the charge Q sector. We have assumed that all the higher energy states are suppressed exponentially at large L t . It is easy to verify that µ ) is the critical chemical potential where the state with charge Q becomes the ground state instead of the state with charge Q − 1. Since these calculations are performed at a fixed lattice size L we define ∆E L (Q) = µ (Q−1) c and these are the numbers reported in Table 3 of the main text.
The computation of µ (Q) c is done by measuring the average charge density as a function of the chemical potential, which is given as c . The use of large L t is essential to suppress the higher excitations such that Eq. (20) provides a reliable fit to the data. We have done further computations to confirm our error analysis. As an example, in Fig 6 we show the data we used to fit to Eq. (20) and extract µ (Q) c and its error from it.

Corrections to scaling
In any scaling analysis of numerical data, corrections to scaling can only be addressed when the leading data does not match the predicted scaling form. In the case when the leading scaling form is obeyed by the data, corrections to the scaling are assumed to be negligible. Otherwise, the fitting becomes an ill-defined problem. Our data for the ratios of the correlation functions, R(L/2) in Eq. (7), are consistent with the leading order scaling form, as shown in the Fig 2 (top and bottom). This suggests that, at the level of accuracy of our data (at the per mille level), the corrections to scaling are negligible.
It is important to note that our results for the conformal dimensions is extracted from much bigger lattice sizes and for a large range of lattice sizes than our results for the Energy. Note that the lattices which are used to extract the conformal dimensions (L/a > 48) are much bigger than L = 8, 10 used in the energy computation. Thus, we do believe that our energy results do contain corrections to scaling. In fact the coefficient c 3 2 = 1.235(10) extracted from the energy computation differs from the coefficient c 3 2 = 1.195(10) by 3%, which we believe is indeed due to scaling violations (or equivalently lattice spacing errors). According to studies in literature [27,28], the leading correction to scaling at the O(2) fixed point is ω ≈ 0.8. Using the relation and substituting c L=∞  Table 2 of Ref. [13]. Following Ref. [13], the ground state at Q = 0 is taken to be zero and the first excited state at Q = 0 is normalized to unity. The open and filled pentagons show results from exact diagonalization studies (ED) for the ground state (GS) and first excited (1st excited) states respectively. Open and filled squares show the same results obtained from epsilon expansion (eps-exp) studies. The circles and triangles are our Monte Carlo (MC) results for the ground state at L = 8 and 10 respectively and are in agreement with previous work. ω = 0.8, we obtain γ ∼ −0.4, which seems like a reasonable estimate of the scaling correction term. Further, also note that in the ratio of correlation functions, the leading order scaling corrections get canceled: For lattices L/a > 48 this correction is less than the statistical errors reported in Table 2 of the main text, and explains why the scaling corrections were not essential in the extraction of conformal dimensions.
Comparison with Results of [13] Recently, there has been some interesting work (see [12,13]), in which the authors compute the energy spectrum of an O(2) symmetric Hamiltonian for different global charge sectors on a finite torus, using both the analytic epsilon-expansion (eps-exp) technique, as well as numerical calculations using exact diagonalization methods (ED). The results are given in Table 2 of Ref. [13]. It is important to note that in the ED method it is difficult to compute the exact ground state energy as a function of Q, and hence the energy at Q = 0 was set to zero. Also, energy units were chosen such that the energy gap to the first excited state in the Q = 0 sector is normalized to unity. With these choices the spectrum computed in [13] in the two methods is plotted in Fig 7. It is very useful to compare our results for the ground state with these results. If we normalize our results similarly, our values for the ground state energies for both L=8,10 lattices, show good agreement with the work of [13] as shown in Fig 7.