Mass anomalous dimension of Adjoint QCD at large N from twisted volume reduction

In this work we consider the $SU(N)$ gauge theory with two Dirac fermions in the adjoint representation, in the limit of large $N$. In this limit the infinite-volume physics of this model can be studied by means of the corresponding twisted reduced model defined on a single site lattice. Making use of this strategy we study the reduced model for various values of $N$ up to 289. By analyzing the eigenvalue distribution of the adjoint Dirac operator we test the conformality of the theory and extract the corresponding mass anomalous dimension.


Introduction
The SU (2) gauge theory with two adjoint Dirac fermions, known as Minimal Walking Technicolor (MWT) [1,2], has been the subject of many lattice studies, all of which have found it to be a conformal theory with a fairly small mass anomalous dimension γ * [3,4,5,6,7]. The most recent measurement obtained by fitting the mode number of the Dirac operator gave a very precise value [8]. The mode number method has also been used to follow the running of γ over a range of energy scales for the SU (3) theory with many light fundamental fermions [9].
The large N version of MWT, the SU (N ) gauge theory with two adjoint fermions, is interesting for several reasons. From a phenomenological point of view, it is expected to be similar to the SU (2) theory. For example, the universal first two perturbative coefficients of the beta function are independent of N . Hence, as in the case of MWT, they point towards the existence of an infrared fixed point with a mass anomalous dimension that is also independent of N . Moreover, numerical results for the mass anomalous dimension for the SU(2) and SU(3) theories appear to be in agreement [10], suggesting that this Nindependence may be a good approximation all the way down to N = 2. From a more theoretical point of view, the large N theory is better suited for connecting with results obtained from different approaches, such as the AdS/CFT correspondence. Fortunately, the numerical study of the infinite volume theory at large N is made possible by the concept of large N volume independence. This implies the equivalence with a single site lattice reduced model, for which simulations can be performed at large values of N , that would be prohibitively expensive on a conventional L 4 lattice. In this context, the study of large N Yang-Mills theory with adjoint fermions has attracted much attention [11]- [20].
In this work we will be using the twisted reduction technique [21,22]. For the adjoint fermion case the specific form of the action has been given in Ref. [23]. This model has been shown to lead to a softer N dependence than the Adjoint Eguchi-Kawai model with periodic boundary conditions [12,23]. The twisted model depends on the choice of the twist tensor. Here we will follow the same symmetric twist prescription as for the pure gauge theory, in which N is taken as the square of an integer number N =L 2 , and the flux through each plane is equal to ±kL (modulo N ). With appropriate values of the integer k, this choice has proven effective in avoiding symmetry breaking for the pure gauge theory [24]. An important advantage of the twisted reduction method is that the dominant 1/N corrections amount to finite size effects on anL 4 lattice. This allows an estimate of the values of N at which the simulations should be performed. In this work we will be using values of N up to 289, corresponding to lattices of size 17 4 .
In summary, the purpose of this paper is to analyze the behaviour of the SU (N ) gauge theory with two flavours of adjoint fermions in the large N limit. Our main goal is to determine whether the theory has indeed a non-trivial infrared fixed point (IRFP) and to measure the mass anomalous dimension at this fixed point. In previous papers some of the present authors studied the behaviour of Wilson loops and the corresponding string tension [25,26,27,28]. Although the results were consistent with the conformal behaviour characteristic of an IRFP, the extraction of the mass anomalous dimension had large uncertainties. Our methodology here will be based on an alternative procedure which has produced very precise estimates in the study of MWT [8]. Preliminary results have been presented in Refs. [29], [30].
The strategy is to determine the anomalous dimension from the structure of the eigenvalue density ρ(ω) of the massless Dirac operator / D. The eigenvalue density is defined as where the sum runs over all eigenvalues iω k of / D. In a mass-deformed conformal field theory (mCFT), this quantity should vanish for w → 0 as [31] where γ * is the mass anomalous dimension at the infrared fixed point, V is the lattice volume, and m is the mass. This behaviour should be contrasted with the one characteristic of a chirally broken theory where the eigenvalue density does not vanish at the origin. In this paper we will use the previous idea to determine γ * for the large N gauge theory with two adjoint quarks. It is clear from the previous formula that it is crucial to work in the region of very small masses and keeping finite volume effects under control. As a reference we will compare our result with those obtained for the pure gauge theory (n f = 0), for which the eigenvalue density has the characteristic behaviour of a chirally broken gauge theory.
The structure of the paper is as follows. In the next section we will collect all the technical aspects concerning the simulation and the extraction of the mass anomalous dimension from the data. In the following we will present the results of our analysis. The paper ends with the presentation of our conclusions.

Methodology
As explained in the introduction, our approach to the large N limit is based on reduction. Hence, we simulate the twisted reduced SU (N ) model on a single site with two adjoint Dirac fermions [23]. In the large N limit the theory is equivalent to the infinite volume lattice gauge theory. For finite N =L 2 , the corrections amount to finite volume corrections in an  Table 1: Number of configurations used for calculating the eigenvalue spectrum at each value of b, κ, N and flux k. The integer parameterk satisfies kk = 1 (mod √ N ). All configurations are separated by at least 25 molecular dynamics updates, and at N = 289 they are separated by 125 molecular dynamics updates.
As explained in the introduction, we have chosen the symmetric twist configuration, with values of the flux integer parameter k given in Tab. 1. These fulfill the condition k/ √ N > 1/9 which was found necessary for the pure gauge theory (TEK model) to respect the center symmetry [24]. This symmetry is a necessary ingredient in the proof of reduction by Eguchi and Kawai. The addition of light adjoint fermions should help in preserving the symmetry but, as shown in Ref. [23], adhering to the condition allows the study of the full range of κ values and leads to a smoother N dependence. As an example of the behaviour of Polyakov loops, which act as order parameters of the center symmetry, in Fig. 1 we display the expectation value of the modulus of the unit winding loop as a function of N .
By definition, this quantity is always positive, but as seen in the figure its size decreases with N for all values of κ.  For each configuration we compute the low-lying spectrum of the modulus square of the massive lattice Wilson Dirac operator in the adjoint representation. This operator is positive definite and in the naive continuum limit corresponds to a 2 (− / D 2 + m 2 ). Thus, its eigenvalues, labelled (aΩ) 2 , are related to those of / D by the expression where a is the lattice spacing and m is the adjoint quark mass. The lowest eigenvalue of our lattice operator defines the spectral gap. In the continuum it is bounded from below by a 2 m 2 . In the case of QCD and for the lattice Wilson Dirac operator, the median of the gap distribution was found empirically [32] to satisfy the relation so that aΩ 0 is proportional to am. In our case, however, we expect that the bound is not saturated at finite N . The reason being the absence of zero-momentum quark states in the reduced model. Hence, quarks are then produced with at least the minimum momentum 2π/L, whose square decreases linearly with 1/N . Since our main goal is the determination of the mass anomalous dimension with the idea presented in the introduction, it is important to study the region close to the critical point and for large values of N to minimize small effective volume corrections. For that reason, our main analysis was based on the study of the lowest 2000 eigenvalues (aΩ) 2 at N = 289 and the lowest 1000 eigenvalues at N = 121. The set of values of b and κ were given in Tab. 1. The computational cost increases considerably as we approach κ c , explaining the smaller number of configurations for that case. Fortunately, the distribution of the lowest lying spectrum does not seem to fluctuate strongly at those values of N .
In determining the value of γ * from the distribution of eigenvalues there are certain alternative procedures which we will describe below.

Determining γ * from a fit to the spectral density
In the continuum γ * could be determined by fitting the spectral density to the form expected for a mass-deformed conformal field theory, Eq. (1.2). However, in order to compare to the lattice data, it is more convenient to look at the spectral density of − / D 2 , given by: This quantity is obtained on the lattice by counting the number of eigenvalues of the modulus square of the Wilson Dirac operator within a bin of size ∆ 2 around (aΩ) 2 . Representing this number by N (aΩ, ∆), the lattice spectral density is given by: (2.4) where N 2 represents the lattice volume on the reduced lattice. The continuum formula for ρ(ω 2 ) gives the following parameterisation for the lattice data: allowing the determination of the three free parameters: B, (am) 2 and γ * . The lowest part of the eigenvalue distribution is the one most affected by finite volume and finite mass effects, hence the fits have to be performed in an intermediate range of eigenvalues aΩ min < aΩ < aΩ max , which preserves the separation of scales on the lattice, From a practical viewpoint the 3-parameter fit demands very precise data in a wide range of eigenvalues and induces strong correlations between the parameters. In some cases the number of parameters can be reduced by assuming that the mass is negligibly small (am = 0). Alternatively we can continue to work at finite mass but use information coming from the smallest eigenvalue to fix the parameter (am) 2 of the fit.

Determining γ * from a fit to the mode number
Alternatively we can follow the procedure introduced in Ref. [8] and extract γ * from the mode number ν(Ω) of the Dirac operator. It is simply defined as the number of eigenvalues of (− / D 2 + m 2 ) below some value Ω 2 . Hence, it is given by V times the integral of the eigenvalue density. We can split this integral into two parts as follows The first part contains the range of eigenvalues which are more sensitive to finite volume and/or finite mass effects. For the second part we can insert Eq. (1.2) and perform the integration to give To determine ν(Ω) on the lattice we just simply count the number of eigenvalues of the modulus square of the lattice Wilson massive Dirac operator below some value (aΩ) 2 . The continuum formula for the mode number implies that the lattice data can be parametrized as follows: where Eq. (2.9), normalised dividing by the lattice volume, is the expression used in Ref. [8] to fit the lattice data. This allows the determination of its four free parameters (ν 0 , A, (am) 2 and γ * ). For the same reason as discussed before for the spectral density, one can attempt to reduce the number of parameters in the fit. If we assume that finite volume (i.e. finite N ) effects are negligible we might take Ω IR close to threshold making ν 0 negligibly small compared to the other term. This leads to a simplified expression which can be used to fit its three parameters (A, am, γ * ) to the modenumber data in a range of eigenvalues Ω ∈ [Ω min , Ω max ] satisfying: If we set am = 0 in Eq. (2.11) we can reduce the number of free parameters even further. This is for example the strategy adopted in Ref. [9] where a 2-parameter fit to the mode number is used, The sensitivity to the fit function, the volume, the mass parameter, or the fitting range will be used to estimate the systematic error in the determination of the mass anomalous dimension.

Analysis of the two flavor case
In this section we will present the results of our study. Our first step is the analysis of the spectral gap of the hermitian Wilson-Dirac operator. To study how this quantity behaves as a function of N and κ, we measure the lowest eigenvalue at b = 0.36 on a range of configurations for N = 25 − 289 and κ = 0.130 − 0.170. As argued in the previous section, in the twisted model we expect the median of the gap distribution (aΩ 0 ) 2 to differ from (am) 2 by a finite volume correction which, interpreted as non-zero momentum contribution, should depend linearly on 1/N . Indeed, a linear fit of this kind seems to describe our nonperturbative data quite well, as shown in Fig. 2. This allows us to determine the mass parameter (am) 2 in the large N limit, listed in Tab. 2.
We have analysed the probability distribution of the spectral gap p(aΩ 0 ) of the Wilson Dirac operator. In the large volume limit of QCD, the analysis in Ref. [32] showed that  the distribution is gaussian: with median proportional to the bare current quark mass and width scaling with the volume and the lattice spacing approximately as ξa = a 2 / √ V . In the reduced lattice this would imply ξa = 1/N . Fig. 3 shows p(aΩ 0 ) for N = 289 at b = 0.35 and 0.36 and several values of κ. A gaussian fit describes the data well. The fitted distribution widths multiplied by N are displayed in Fig. 4 for our n f = 2 configurations at several values of N and κ. Our results follow rather well the behaviour also reported in QCD.
Let us now move on to describe our results for the distribution of eigenvalues. One of the main points is to analyze the N -dependence of this distribution. We already saw that this dependence affects the gap of the spectrum, but we expect this effect to have a small impact for higher eigenvalues. This can be seen in Fig. 5, which shows a histogram of the number of eigenvalues as a function of (aΩ) 2     We have seen that, at least for the two biggest values of N , the eigenvalue distribution roughly coincides beyond a certain threshold value. The question now is to see if this distribution behaves as expected from the IRFP hypothesis and to extract γ * from it. In the previous section we gave two alternative methods of fitting the data. One is to compare the eigenvalue distribution with Eq. (2.5). The other is to compute the mode number and fit it to Eq. (2.8) or to its simplified expression Eq. (2.11).
In performing a fit one has to select the range of values Ω min < Ω < Ω max to be fitted. A lower value of Ω min increases the sensitivity to the value of the mass (am) 2 but also risks to be more affected by finite effective volume (finite N ) corrections. For the mode number fit the same is expected to happen for the parameter ν 0 . Furthermore, the narrower the fitting range the stronger the correlations among parameters leading to high uncertainties in γ * .
If we use only the data which is least affected by finite volume and finite mass effects, we can produce our most precise determination of γ * : Hence, we take our data of N = 289 and κ = 0.17 and fit the distribution of eigenvalues to Eq. (2.5) with (am) 2 set to zero. The upper edge of the fitting range (aΩ max ) 2 covers almost all of the N = 289 data. For the lower edge (aΩ min ) 2 the cut is set to twice the lowest eigenvalue for N = 121. The result for b = 0.36 is γ * = 0.268(2) and for b = 0.35 is γ * = 0.271 (1). It is remarkable that both values of b give consistent results within the purely statistical 1% errors. To give an idea of the quality of the fit we display the data in Fig. 6 together with the best fit. In the figure we also include the data of N = 121 at both values of b. The fitted function also describes well the behaviour of the data at N = 121, except for the smallest eigenvalues where finite volume effects should be mostly felt. The continuous line going through the κ = 0.165 data was obtained fitting the corresponding data with γ * fixed to the value obtained at κ = 0.17 and the mass square (am) 2 to the large N extrapolated value given in Tab. 2.
The analysis of the previous figure and fits indicates that our data look consistent with the predictions of an infrared fixed point with mass anomalous dimension close to 0. 27 Figure 7: Fitted γ * as a function of (aΩ min ) 2 from fitting the spectral density data to Eq. (2.5). For each point, a range of masses (am) 2 between zero and the lowest eigenvalue at N = 289 are used, and the upper fit range bound (aΩ max ) 2 is also varied over a range of values from 0.25 up to the largest value for which we have data. The black line shows our best estimate for γ * with statistical error only.
uncertainties involved in our analysis and to compare our results with those corresponding to theories which are not conformal in the IR. The latter will be done in the next section where we will compare our n f = 2 results with those obtained for n f = 0. What we will now present is an evaluation of the systematic errors. This point will be analysed by estimating the effect of finite N corrections, finite mass corrections and sensitivity to the fitting range. We will also test what results come out if we use the mode number instead of the eigenvalue density. A good summary of the effect of all systematics on the value of γ * is provided by Fig. 7. Here we display different determinations of γ * for the N = 289 and N = 121 data by varying Ω max , Ω min and (am) 2 . The latter is varied within the range extending from zero to the minimum N = 289 eigenvalue. The data is displayed as a function of (aΩ min ) 2 . For a large range of x-axis values all determinations of γ * (obtained by using different values of Ω max and (am) 2 ) fall within a horizontal strip whose width serves as an upper bound to our systematic error δγ * = 0.05. For larger values of (aΩ min ) 2 the fitting range narrows, obviously leading to a wider spread of values of γ * .
To better analyze the dependence of the fitted value of γ * with the remaining parameters, we display in Fig. 8 its correlation with the mass parameter used in the fits for all values (aΩ max ) 2 and (aΩ min ) 2 < 0.07. This shows that the main source of systematic errors is indeed the value of the mass.
A similar analysis can be done for the mode number. Results are essentially compatible.     Figure 10: Eigenvalue density distribution, in log-log scale, at b = 0.36 for n f = 2 and n f = 0. The lines correspond to the fits to Eq. (2.5) described in the text.
example Fig. 9 shows a fit to Eq. (2.11) using the same configurations and fit ranges as for the eigenvalue density in Fig. 6, which gives γ * = 0.270(2) for b = 0.36 and γ * = 0.272(1) for b = 0.35, which are in good agreement with the numbers determined from the eigenvalue density.

Comparison between zero and two flavors
The comparison between our n f = 2 results with those obtained for n f = 0 is essential to substantiate the claim of conformality. The latter theory is not conformal and should display a different behaviour. Indeed, for the pure Yang-Mills theory it makes no sense to speak of γ * itself, since there is no infrared fixed point. However, it still makes sense to study the behaviour of the spectral density and mode number distribution as a function of its argument.
To carry on the previous study, we generated pure gauge configurations at the same values of b and N and a range of κ values, and repeated our previous analysis for this data.
The κ values were chosen in such a way as to explore similar small values of the minimum eigenvalue.
The fits to the spectral density for the largest values of κ were qualitatively as good as those of n f = 2. Fig. 10 shows, in log-log scale, the fits to the lattice spectral densities for n f = 2 and n f = 0 at b = 0.36. The values of γ * in this plot are obtained by fitting the spectral densities to Eq. (2.5). The mass (am) 2 is varied in the range 0.5(aΩ 0 ) 2 ≤ (am) 2 ≤ (aΩ 0 ) 2 , the number of bins is varied from 40 to 80, and the lower edge of the fit range is varied in the range 1.5(aΩ 0 ) 2 ≤ aΩ min ≤ 2(aΩ 0 ) 2 . The final point is an average over all of these choices for many bootstrap replicas of the data. Keeping γ * fixed to this average, the fits are repeated for each set of data in order to determine am. For this final fit the lower edge of the fit range is set to 2(aΩ 0 ) 2 . In all cases the χ 2 per degree of freedom of the fit is below 2. Hence, from the quality of the fit to the spectral density one cannot deduce the presence of an infrared fixed point.
One remarkable difference appears when looking at the dependence of the extracted value of γ * with the mass (am) 2 . The result is shown in Fig. 11a. The value of γ * seems quite stable for the n f = 2 case. This is what one expects in the vicinity of an infrared fixed point as the anomalous dimensions tends to a constant at the fixed point. The result for n f = 0 is quite different, showing a pronounced drop as we move away from the critical value of κ. For the smallest masses the fitted value of γ * reaches as high values as 0.8. Notice that a value equal to 3 would imply a constant value of the spectral density at the origin. Our data show a growing γ * for lighter masses but do not yet reach the value of 3 predicted by chiral symmetry breaking, through the Banks-Casher formula. This is probably due to finite volume and/or finite mass effects.
Another marked difference between the two cases appears in the dependence of the spectral densities on the bare inverse coupling b. In the conformal case one expects the coupling to be a marginally irrelevant operator close to the IRFP, in contrast to the n f = 0 case where it is marginally relevant and determines the lattice spacing. The values of γ * displayed in Fig. 11a for n f = 2 do indeed have a very small dependence on b. However, this dependence is large for the case of n f = 0. In the quenched case, the results corresponding to the two different couplings only show scaling when expressed in physical units in terms of the string tension. This is shown in Fig. 11b, where we have used the data for the n f = 0 string tension obtained in Refs. [33] (σa 2 = 0.09 and σa 2 = 0.043 for b = 0.35 and b = 0.36 respectively).

Conclusions
We have performed a measurement of the mass anomalous dimension γ * of the SU(N ) gauge theory with two adjoint Dirac fermions, in the large N limit using the concept of large N twisted reduction. Results from a single site lattice model at large values of N , have the expected qualitative behaviour of the spectral density and mode number of the adjoint massless Dirac operator. The distribution for small masses (extracted from the lowest eigenvalue) can be well-fitted with the expectations of an IRFP. From the data we extract a value of γ * = 0.269±0.002±0.05, where the first error is statistical and the second one systematic. This value is similar to previous lattice determinations of this quantity for the SU(2) theory.
Does our result provide conclusive evidence of the presence of an infrared fixed point for the SU(∞) gauge theory with 2 flavours of adjoint fermions? To try to answer this question we repeated the analysis for the n f = 0 case, which is known to have a completely different behaviour at criticality. However, we observe that the spectral densities at fixed b and small quark mass, can also be fitted with the same formulas with a larger value of γ * . This conveys a word of warning about drawing conclusions about the existence of an IRFP only from the capacity to fit the spectral density or mode number with a powerlike distribution. Nevertheless, there are marked differences between the behaviour observed in the n f = 0 and n f = 2 cases. One of them is the dependence of the result on the bare coupling b. The n f = 2 results for our two values of b, 0.35 and 0.36, are consistent with each other. This is not the case for the n f = 0 data. Although insensitivity to the bare coupling b is certainly the expected result for an IRFP, it is difficult to exclude the fact that this is not simply due to the smaller value of the beta function when adding fermions in the adjoint. The second difference refers to the change of behaviour as we approach criticality κ −→ κ c . The extracted value of γ * for n f = 2 remains fairly stable as one expects if the behaviour is indeed dictated by the presence of an IRFP. On the contrary for the n f = 0 data we observe a pronounced rise of the value of γ * as we approach criticality.
To improve on these results using this method, smaller fermion masses and larger volumes (i.e. larger values of N ) would be required. Since finite volume effects are the dominant limitation, one possibility for future work would be to extend the single site lattice to a 2 4 lattice.