Conformal dimensions in the large charge sectors at the O(4) Wilson-Fisher fixed point

We study the O(4) Wilson-Fisher fixed point in 2+1 dimensions in fixed large-charge sectors identified by products of two spin-j representations $(j_L, j_R)$. Using effective field theory we derive a formula for the conformal dimensions $D(j_L, j_R)$ of the leading operator in terms of two constants, $c_{3/2}$ and $c_{1/2}$, when the sum $j_L + j_R$ is much larger than the difference $|j_L-j_R|$. We compute $D(j_L,j_R)$ when $j_L = j_R$ with Monte Carlo calculations in a discrete formulation of the O(4) lattice field theory, and show excellent agreement with the predicted formula and estimate $c_{3/2}=1.068(4)$ and $c_{1/2}=0.083(3)$.


INTRODUCTION
Conformal field theory (cft) holds a central place in the study of quantum field theory (qft), as it is relevant to both particle physics and condensed matter systems at criticality, and via the gauge/gravity correspondence even to the description of quantum gravity. Generically, cfts do not contain any small couplings that can be used in a perturbative analysis. However, the conformal symmetry constrains its observables such that we can determine any n-point function using only operator dimensions and threepoint function coefficients. While it is possible to treat strongly coupled theories with methods such as the large-N expansion, the -expansion (see [1] for a review) and the conformal bootstrap [2], they are notoriously difficult to access analytically. In simple cases, Monte Carlo (mc) techniques offer a reliable numerical alternative [3,4].
Recently, it has been shown in a series of papers [5][6][7][8][9][10] that working in a sector of large global charge results in important simplifications and gives us a perturbative handle to study cfts using effective field theories (efts): it is possible to write an effective action as an expansion in terms of a large conserved charge with unknown coefficients. For the Wilson-Fisher point in the three-dimensional O(N ) vector model [11], except for two low-energy couplings, all terms are suppressed by inverse powers of the large charge [6]. The approximate physics of the cft becomes accessible as a function of these two couplings which we label as c 3/2 and c 1/2 . This suggests a double-pronged approach to cfts, which involves using the large-charge expansion to determine the effective action, paired with mc calculations to determine the low-energy couplings. For the case of the O(2) Wilson-Fisher cft, this approach has been successfully implemented recently [12]. In particular, it was shown that the predictions obtained with the two couplings remain very accurate even for low charges.
In this letter, we explore the viability of this approach for the O(4) Wilson-Fisher cft, which has qualitatively distinct features from the O(2) model studied earlier. The fact that O(4) symmetry is non-Abelian and that it leads to two conserved global charges j l and j r , creates novel challenges. The ground state can become spatially inhomogeneous requiring a different analysis in the eft, and the construction of a worldline-based lattice model becomes necessary to access easily the large-charge sectors. The cft with O(4) symmetry is also interesting in many subfields of physics. For example, it arises naturally in the study of finite-temperature chiral phase transitions in two-flavor quantum chromodynamics (qcd) with massless quarks [13,14]. It is also of interest in studies of strongly correlated electronic systems at half filling built out of models of interacting electrons with spin [15].
Traditional O(4) lattice models are constructed using classical vectors. Unfortunately, in the study of large charge sectors, using traditional mc methods based on sampling classical vectors leads to severe signal-to-noise ratio problems. While worldline representations can in principle solve these problems [16,17], the presence of an infinite Hilbert space at each lattice site can still lead to algorithmic inefficiencies. Fortunately, a discrete version of the O(4) model with a finite Hilbert space per lattice site is easy to construct [18,19]. In this work we use this discrete formulation to accurately compute the conformal dimensions D(j l , j r ) (defined below) at the O(4) Wilson-Fisher fixed point, when j l = j r = j (see Figure 1). The large-charge prediction (see Eq. (5)) is an excellent fit to the lattice data even up to the smallest charge giving us c 3/2 = 1.068(4) and c 1/2 = 0.083 (3).

LARGE CHARGE PREDICTIONS
The eft approach to the O(4) cft is based on the construction of an effective action at large charge. Using the fact that SU l (2) × SU r (2) is a double cover of O(4), this action can be put into the form [6,8] where g(r, t) ∈ SU (2), r is the coordinate on Σ, dg 2 = Tr ∂ µ g † ∂ µ g , and c 3/2 and c 1/2 are the two leading lowenergy couplings referred to earlier. This action is to be understood as an expansion around the fixed-charge ground state. We study the system on a spatial Riemann surface Σ with scalar curvature R. The field g transforms as g → V l gV −1 r under the global SU (2) l × SU (2) r . The corresponding Noether charges are the two-by-two Hermitian traceless matrices in the su(2) algebra: √ 2c 3/2 dg . Their two eigenvalues are ±j l and ±j r . Under the action of SU (2) l × SU (2) r the charges in Eq. (2) transform as but j l and j r remain invariant. We will refer to the class of configurations g connected by SU (2) l × SU (2) r transformations as the (j l , j r ) sector. In the underlying qft the sectors (j l , j r ) naturally label the irreducible representation space of SU (2) l ×SU (2) r , hence the values of j l and j r are quantized, i.e., j l , j r ∈ 1 2 Z. Their sum must be integer, j l + j r ∈ Z, because we consider only states that are representations of O(4).
Instead of Q L,R it is more convenient to work with the projections q l,r = Tr(Q l,r σ 3 )/2. In a fixed (j l , j r ) sector these projections will take values in the range −j l,r ≤ q l,r ≤ j l,r . It is natural to identify them with the quantized charges of the states in the representation with highest weights (j l , j r ). In Ref. [6] it was shown that the minimal energy solutions to the equation of motion (eom) for the action (1) for fixed values of q l,r are homogeneous in space and arise in sectors with j l = j r = max(q l , q r ). This leads to the formula for the minimal energy in a fixed (j, j) sector: where ζ(s|Σ) is the ζ-function for the Laplacian on the surface Σ and represents the contribution of the Casimir energy. Since in a cft the conformal dimension of an operator is identified with the energy on the unit sphere Σ = S 2 of the corresponding state, we deduce a formula for the dimension of the lowest operator in the (j, j) sector: where c 0 = ζ(−1/2|S 2 ) ≈ −0.094 is a universal constant [20,21]. In this work we generalize this formula to any representation (j l , j r ). We need to find the minimal-energy solutions admitted by the action Eq. (1) whose charge matrices Q l and Q r are diagonal and correspond, by the argument above, to highest-weight states in a representation of SO(4). In order to study the sectors with j l = j r , the analysis in Refs. [6,22,23] suggests that we need to look for inhomogeneous field configurations g(r, t) of the form g(r, t) = cos(p(r))e iµ1t sin(p(r))e iµ2t − sin(p(r))e −iµ2t cos(p(r))e −iµ1t , where µ 1 and µ 2 are constants parameterizing the action of SU (2) l × SU (2) r on the inhomogeneous configuration encoded by the undetermined function p(r). The homogeneous solution of Ref. [6] with µ 1 = µ 2 = O j 1/2 and p(r) = 0 describes the j l = j r sectors. In order to explore solutions with small non-zero values of |j l − j r |/ max(j l , j r ) we may expand the action in a series in η 2 = (µ 2 − µ 1 )/(µ 2 + µ 1 ) which will naturally be small. At leading order in η 2 the eom is an elliptic sine-Gordon equation where is the Laplacian on Σ and Λ 2 = µ 2 2 − µ 2 1 . The case Σ = T 2 was already discussed in [22,23]. Here we concentrate on Σ = S 2 (r 0 ) in order to calculate operator dimensions. We express the parameters µ 1 and µ 2 as functions of the eigenvalues j l and j r and find that the leading (tree-level) contribution to the energy of the solutions to the eom has the same form as Eq. (4), where now j = j m = max(j l , j r ), plus an extra contribution that captures the inhomogeneity of the solution (i.e., ∇p(r) = 0): As discussed in the supplementary material, when Σ = S 2 (r 0 ), the eom (7) admits different branches of smooth solutions, parameterized by an integer which counts the zeros of p(r). The energy is minimal in the first non-trivial branch ( = 1), where 2 ≤ r 2 0 Λ 2 < 6. Here the integral of the divergence can be computed numerically in terms of an expansion in |j l − j r |/j m to give with λ 2 ≈ 0.2455. This is the leading contribution in the large-charge expansion. There will be in general higherorder corrections suppressed by inverse powers of the large charges due to sub-leading terms in the tree-level action in Eq. (1) and to quantum corrections.
There is only one term of order O j 0 : the Casimir energy of the Goldstones resulting from the spontaneous symmetry breaking SO(3) × D × SO(2) 2 → SO(2) × D discussed in the supplementary material. The two broken generators of the isometries on the sphere only give rise to one Goldstone degree of freedom (dof). Together with the 2 dof from the broken internal symmetries, they are arranged into one type-I and one type-II Goldstone field in the notation of [24]. Only the former contributes to the Casimir energy as . The zero-point energy is different from the one in the (j, j) sector because the low-energy excitations only propagate in the direction of the unbroken sphere isometry. Once again we can use the state/operator correspondence and obtain the final formula for the conformal dimension of the lowest operator in the representation (j l , j r ) of SO(4) when j l = j r : As we have stressed, the conformal dimensions only depend on the two Wilsonian couplings c 3/2 and c 1/2 , which are the same coefficients that appear in Eq. (5) for the j l = j r case. We now explain how we determine them using mc methods with our lattice model.

LATTICE SIMULATIONS
Our lattice model was first introduced in Ref. [18] as a model for pion physics in two-flavor qcd and studied with an efficient mc algorithm. It is constructed using four Grassmann fields ψ α (x), ψ α (x), α = 1, 2 at every three-dimensional periodic cubic lattice site x = (r, t) of size L in all the directions. If we arrange these four-fields into a 2 × 2 matrix of the form g αβ (x) = ψ α ψ β we can write the lattice action as where xy are nearest-neighbor bonds. This action is on even sites. The partition function of the model can be expressed as a sum over configurations where each site either contains a vacuum site or a worldline of an O(4) particle in the vector representation. Thus, each worldline has four possible states that label the eigenvalues (q l , q r ) = (±1/2, ±1/2) of particles that travel through the sites. These can be thought of as oriented loops with two colors (say red and green). An illustration of a configuration is shown in Fig. 2.
The weight of a worldline configuration is given by U Nm where N m is the number of vacuum sites. As U is tuned, the model undergoes a phase transition between the massive symmetric phase at large values to a phase where the O(4) symmetry is spontaneously broken at small values. Using well-established mc methods [18,19] we first demonstrate that at the critical point we obtain the O(4) Wilson-Fisher cft by computing the critical exponents ν and η. For this purpose we compute the current susceptibility ρ s and the order parameter susceptibility χ, details of which can be found in the Supplementary Material. Finite-size scaling theory near a second-order phase transition predicts that ρ s L and χL 2−η must be simple polynomials of (U − U c )L 1/ν . A combined fit of our data gives U c = 1.655394(3), ν = 0.746(3) and η = 0.0353(10). In Fig. 3 we plot our data and the fit. These exponents are in excellent agreement with earlier results, ν = 0.749(2) and η = 0.0365(10), obtained from the traditional lattice model [25].
Having established that our lattice model indeed reproduces the O(4) cft when U = U c , we can use the method we developed in Ref. [12] to accurately compute the conformal dimensions D(j, j) at the O(4) cft. We can create configurations in a specific (j l , j r ) sector by placing appropriately charged sources and sinks at t = 0 and t = L/2 respectively. More concretely, sources that create a red loop are assigned the charge (1/2, 1/2) and the sinks that annihilate them are assigned the charge (−1/2, −1/2). Similarly, those that create and annihilate the green loops are assigned charges (1/2, −1/2) and (−1/2, 1/2). Using these fundamental sources we can construct sources and sinks with any charge (q l , q r ). However, since each site can only have one red or one green source, to create a source with a large charge we distribute the fundamental sources in a local region near the origin (see Supplementary Material for more details). Since the couplings c 3/2 and c 1/2 can be computed by fitting the data j D(j, j) j D(j, j) (this work) (from [26] I. Results for the conformal dimensions D(j, j) up to j = 5 computed using worldline mc methods in this work (second and fifth column). We also compare our results with earlier calculations up to j = 2 found in [26].
for D(j, j) to the predicted form in (5), in this work we only study the sector with j l = j r = j, j = 1/2, 1, 3/2, . . . . For this purpose we only work with sources and sinks of equal charges by creating 2j sources of red loops at t = 0 and annihilating them at t = L/2. This naturally projects us into the highest-weight representation sector with j l = j r = j. Let Z j (L) be the partition function in the presence of these sources and sinks. In Ref. [12] we developed an efficient algorithm to compute the ratio R j (L) = Z j (L)/Z j−1/2 (L), which is expected to scale as C/L 2∆(j) for large values of L. By evaluating R j (L) for various values of j, L and fitting to the expected form we can accurately compute the difference in the conformal dimensions ∆(j) = D(j, j) − D(j − 1/2, j − 1/2). From these differences we can also estimate D(j, j), since conformal invariance fixes D(0, 0) = 0. Our final results are tabulated in Table I up to j = 5. As the table shows, our results are also in good agreement with earlier calculations up to j = 2 [26]. Fitting the data in Table I to the large j form in Eq. (5) we obtain c 3/2 = 1.068(4) and c 1/2 = 0.083(3) (see Fig. 1).

CONCLUSIONS
In this letter we provide a new prediction for the anomalous dimensions D(j l , j r ) (see (10)) at the O(4) Wilson-Fisher fixed point in terms of the two couplings that appear in the fixed large-charge effective action (1). Our prediction is valid in the limit of large (j l , j r ) and small |j l − j r |/ max(j l , j r ). We then use a discrete lattice O(4) model to compute the two couplings by fitting the data for D(j, j) to the prediction in Eq. (5) obtained from an earlier work. We also demonstrate that this prediction provides an excellent approximation even at small values of j (see Fig. 1). Our estimate c 3/2 = 1.068(4) and c 1/2 = 0.083(3) can be used in (10) to predict D(j l , j r ) even for j l = j r . While our lattice model can in principle be used to check the validity of these predictions, our method is likely to suffer from signal to noise ratio problems when j l and j r are sufficiently large and different. Discrete lattice models like ours can in principle also be designed for other non-Abelian symmetry groups, thus allowing us to explore the robustness of the large-charge eft method for general cfts. Such extensions are likely to bring new challenges providing a fertile ground for further research.

SUPPLEMENTARY MATERIAL Linearized EOM and symmetry breaking pattern
To our knowledge, there is no known solution to the eom in Eq. (7) in terms of elementary functions. Nonetheless, we can still understand some of its qualitative properties via the study of the associated linear problem, obtained by setting p(r) = 2 (r) and taking the limit → 0. At leading order, satisfies the Laplace equation: (r) + Λ 2 (r) = 0.
As it is well known, when Σ = S 2 (r 0 ) this equation admits smooth solutions only for Λ 2 r 2 0 = ( + 1), = 0, 1, . . .. These are the spherical harmonics (θ, φ) = Y m (θ, φ). Imposing reality for selects m = 0. Moreover, we know that energy minimization requires Λ 2 to take the smallest non-vanishing value Λ 2 = 2/r 2 0 . The final solution in a convenient normalization is then In the non-linear equation (7), the values of Λ are not quantized, but for each value of there is a continuous branch of solutions with ( + 1) ≤ Λ 2 r 2 0 < ( + 1)( + 2). It is convenient to express 2 as function of the charges. At leading order in η, which shows that the linear limit 2 → 0 is the same as requiring |j l − j r | max(j l , j r ). The leading contribution to the integral appearing in the energy formula (9) is then: This explicit solution lets us identify the symmetrybreaking pattern associated to fixing the two representations independently. The minimal-energy configuration is not homogeneous, but a function of the polar angle θ. This breaks the SO(3) isometry of the sphere spontaneously to SO(2). The complete breaking pattern, including the internal symmetries is then with D = D − µ 1 SO(2) 1 − µ 2 SO(2) 2 (see also [7]). The conformal group is broken to the sphere isometries SO(3) times the time translations D, and the global O(4) is broken to SO(2) 2 by fixing the charges. This can be seen as an explicit breaking and does not give low-energy Goldstone dofs. The second breaking is spontaneous, and the two broken generators T a (θ, φ) of SO(3) produce only a single low-energy excitation. This is because the equation [27] c a (φ)T a (θ, φ) (θ) = 0 (17) admits one non-trivial solution for the functions c a (φ). All in all, we have three Goldstone dof. A study of the perturbation theory around (θ) shows one field with linear dispersion relation ω = c s p (type-I Goldstone) and one with quadratic dispersion ω ∝ p 2 (type-II Goldstone). Together, they account for all the dof. The speed of sound for the type-I Goldstone is fixed by three-dimensional scale invariance to be c s = 1/ √ 2. The low-energy excitations only propagate in the direction of the unbroken translation φ, so the corresponding Casimir energy is

Lattice Observables
In order to establish that the discrete model described by the action in Eq. (11) flows to the Wilson-Fisher O(4) cft, we have computed two observables: the current and the order parameter susceptibilities.
The current susceptibility is defined as where J µ (x) is the O(4) conserved current at the site x in any fixed direction µ. In a worldline configuration, such as the one illustrated in Fig. 2, one can choose a surface perpendicular to µ and compute the current flowing out of that surface by counting the directed arrows of a given color with appropriate sign. Red loops or green loops will give the same result due to the O(4) symmetry. The order parameter susceptibility is defined as where a † y creates an O(4) particle and a x destroys it. In terms of the worldline configurations, such creation and annihilation events are introduced as sources and sinks (see Fig 4) and can be sampled during the update of the worm (or directed loop algorithm) as described below. Hence, we can compute χ during such an update.

Worm Algorithms
Worm (or directed loop) algorithms are by now well established and easy to construct. Depending on the way the detailed balance is implemented, they can be constructed in different ways. Some algorithms can be found in the earlier work mentioned in the main paper. In this work we have used two different algorithms to make sure that each one is free of errors. When we compute the physical observables ρ s and χ in order to establish the O(4) Wilson-Fisher fixed point, we used an algorithm which we refer to as ALGO1. This algorithm works in the absence of external charges (sources and sinks). We used a different algorithm, which we call ALGO2, when updating configurations in the presence of sources and sinks. Below we give some details of ALGO2. In fact, in the absence of the sources and sinks, ALGO2 can also be used to compute χ and ρ s , and we have verified that they give the same numerical results as ALGO1 within errors.
In the absence of the sources and sinks, ALGO2 begins in a configuration that contributes to the partition function but samples "worm sectors" with one additional source (tail) at y and sink (head) at x of a given color as illustrated in Fig. 4. Each configuration in the worm sector is given a unique weight W (x, y) depending on the location of the source and the sink. At the end of each worm update the source and the sink disappear and the configuration returns to the partition function sector. The configurations generated during the worm update with a source and a sink help measure the order parameter susceptibility Eq. (20). A brief description of the worm update is given below: Begin: Choose a random initial site y and a random color c (red or green).
1. If y is a vacuum site, then propose to create the configuration in the worm sector with color c with a probability W (y, y)/U . If the proposal is accepted based on a Metropolis accept/reject decision, create a configuration in the worm sector with both the head and the tail located on the site y.
2. If y contains a loop of color c, and this loop enters the site y from x, break the loop between y and x and create a worm configuration with the tail at y and the head at x with a probability W (x, y).
loops contribute to χ when the weights W (x, y) are appropriately divided out. On the other hand configurations in the partition function sector contribute to the observable ρ s L. Additional updates, such as changing the color of loops, as well as their orientation were used to achieve faster decorrelation. The above algorithm can also be used to update worldline configurations in a fixed-charge sector as long as we make sure that the worm updates do not change the location of the sources and sinks. This is easily accomplished if we make sure that the initial site does not contain sources or sinks and if a worm move tries to change their locations it is not accepted. In this work we focus on sectors with charges (q l = q r = j). Since in our lattice model we cannot create all the charges on a single site, we spread them around on a two-dimensional plane around a central site on the t = 0 and t = L/2 slices. One example of how we spread out the charges in the plane on both time slices is shown in Fig. 5. The central site is labeled as 1. For example, a charge with j = 3/2 could be created with sources and sinks of red loops at sites labeled as 2, 3, 4. Similarly, the charge with j = 3 could be created with sources and sinks of red loops at sites 2, 3, 4, 5, 6, 7. We always leave out the central site for a reason that we explain below. We can define the sum over weights of such configurations with sources and sinks separated by L/2 in the sector with charge j as the partition function Z j (L). We can then use our worm algorithm to sample the configurations that contribute to Z j (L). We can then also measure ratios of the form For this purpose we first generate an equilibrium configuration that contributes to Z j−1/2 (L). We then peform a special worm update called a measurement update that creates an additional j = 1/2 charge (red loop) at the central site on the t = 0 sheet. If this loop touches the central site on the t = L/2 we get a contribution to R j (L/2). In order to increase the efficiency of the measurement process we enhance the worm weights W (x, y) ∼ |t xy | p so that the worms do grow large enough to contribute to the observable. This tuning is important and helpful when j becomes large since there the ratio R j (L) is small.
We have checked that the algorithm can efficiently sample large j sectors. We also find that the charges between the source and the sink flow either through the bulk, or over the boundary with equal probability. In addition, we have also checked that the choice of distribution of sources and sinks does not change the final results.

Extraction of conformal dimensions D(j, j)
In the main text, we have briefly stated how the idea initially developed in Ref. [12] can be used to calculate the conformal dimensions D(j, j) for the fermionic realization of the O(4) model. Here we provide further details. In the section above we explained how the measurement update can be used to compute the ratio R j (L). This observable gives us directly access to the conformal dimensions through the relation where Z j (L) is the partition function in the charged sector with q l = q r = j. At the critical point we expect Z j (L) ∼ 1/L 2D(j,j) which implies that R j (L) behaves as C/L 2∆(j) for large lattices, where ∆(j) = D(j, j) − D(j − 1/2, j − 1/2). The term γ L ω are the scaling corrections, which are too small to be detected in our data for the lattices we have used. Hence, a single power law fit was adequate for the extraction of ∆(j, j). The results for R j (L) for various values of j and L are shown in Fig. 6.
The absolute conformal dimensions D(j, j) are obtained using the relation = ∆(j, j) + ∆(j − 1/2, j − 1/2) + · · · + ∆(1/2, 1/2) (23) and using that D(0, 0) = 0. Moreover, since ∆(j) are extracted from fits to independent sets of simulations, there is no correlation between the different values of ∆(j). Hence, the errors in ∆(j, j) can be added in quadrature to obtain the error on D(j, j). The straight line fit on a log-log plot is indicative of the power law behavior, and the slope gives the difference of the conformal dimensions 2∆(j). A very good χ 2 /DOF 1 is obtained with a single power law for data points, indicating negligible scaling corrections, as well accurate extraction of the conformal dimension of the lowest operator in the (j, j) sector. Due to the efficient worm algorithm, there is no visible signal-to-noise problem in these correlators.
Using the extracted values of ∆(j) and D(j, j), we perform a fit to the predicted conformal dimensions in the sector j l = j r = j as given by formula 5. For this, we simultaneously fit the quantities ∆(j) and D(j, j), keeping the value of c 0 = −0.09372 fixed. This yields the values of the coefficients to be c 3/2 = 1.068(4) and c 1/2 = 0.083 (3). The systematic errors, as determined by changing the fit range exceeded the statistical error, and is the main source of the quoted errors. Although c 0 is known analytically, in order to determine how well our data is able to estimate this constraint, we performed a different fit where we allowed all the three coefficients c 3/2 , c 1/2 and c 0 to vary. The values of c 3/2 and c 1/2 obtained with this strategy are consistent with the previously quoted values, but the fit gives c 0 = −0.072 (10), which agrees with the theoretical estimate within two sigmas. The χ 2 /DOF in all the fits was around 0.01 − 0.20, indicating the very good quality of the fits. The presence of the subleading terms in the expansion in Eq. 5 is not observable in our calculations.