Non-perturbative determination of the $\Lambda$-parameter in the pure SU(3) gauge theory from the twisted gradient flow coupling

We evaluate the $\Lambda$-parameter in the $\overline{\mathrm{MS}}$ scheme for the pure SU(3) gauge theory with the twisted gradient flow (TGF) method. A running coupling constant $g_{\mathrm{TGF}}^2(1/L)$ is defined in a finite volume box with size of $L^4$ with the twisted boundary condition. This defines the TGF scheme. Using the step scaling method for the TGF coupling with lattice simulations, we can evaluate the $\Lambda$-parameter non-perturbatively in the TGF scheme. In this paper we determine the dimensionless ratios, $\Lambda_{\mathrm{TGF}}/\sqrt{\sigma}$ and $r_{0}\Lambda_{\mathrm{TGF}}$ together with the $\Lambda$-parameter ratio $\Lambda_{\mathrm{SF}}/\Lambda_{\mathrm{TGF}}$ on the lattices numerically. Combined with the known ratio $\Lambda_{\overline{\mathrm{MS}}}/\Lambda_{\mathrm{SF}}$, we obtain $\Lambda_{\overline{\mathrm{MS}}}/\sqrt{\sigma} = 0.517(10)(^{+8}_{-7})$ and $r_{0}\Lambda_{\overline{\mathrm{MS}}}=0.593(12)(^{+12}_{-9})$, where the first error is statistical one and the second is our estimate of systematic uncertainty.


Introduction
The Λ-parameter is a fundamental quantity in asymptotically free gauge theories and plays the role to set the scale of the theory. Λ characterizes the low energy non-perturbative physics and its determination is one of the most important tasks in lattice gauge theory. In the pure Yang-Mills theory, Λ is the only free parameter of the theory and is determined from the coupling constant. Its value depends on the renormalization scheme. In the MS scheme, for example, it is defined by , for the pure SU(N C ) gauge theory. Since the MS scheme is only defined perturbatively, the non-perturbative estimate of Λ MS thoroughly within the MS scheme is impossible. Therefore we usually convert a Λ-parameter determined with a non-perturbative scheme to Λ MS through the perturbative relation.
On the lattice, the Λ-parameter can be defined by with the lattice spacing a. The bare coupling g 0 can be related to the lattice spacing a non-perturbatively and be used as g Lat (1/a) = g 0 in eq. (1.2). This defines a lattice scheme. It is, however, well known that the scaling is largely violated for the range of g 2 0 accessible with the presently available computational power. In the early stage of the lattice studies, it was common to use an improved coupling such as g 2 E = 8N C N 2 C −1 (1 − u p ) or g 2 A = g 2 0 /u p with u p the observed plaquette value [1,2]. They exhibit a better scaling property, nonetheless there are only intuitive arguments of "tad-pole improvement" to explain why they work.
Great progresses for evaluating non-perturbatively running coupling constants have been made with the discovery of the step scaling method [3], where the renormalization scale is introduced by the physical box-size of the target system. In this method, one can calculate the running coupling in a wide range of the scale covering both the hadronic scale, where we make the non-perturbative calculation of physical quantities with lattice techniques, and the high energy scale, where we can estimate the Λ-parameter neglecting higher order corrections. The most successful non-perturbative scheme for the running coupling constant in QCD is the Schrödinger functional (SF) scheme [4][5][6][7][8][9], in which a specific Dirichlet boundary condition is imposed on the temporal direction of the box. The advantages of the SF scheme are that it is regularization independent and can be defined non-perturbatively. In addition, the calculation of the Λ-parameter ratio Λ MS /Λ SF has been done in ref. [10] perturbatively. The disadvantage of the SF scheme, on the other hand, is that it becomes difficult to calculate the coupling at larger physical box sizes (i.e. low energy renormalization scale) due to the appearance of exceptional configurations and the noisy behavior which result the large statistical error [6].
Several other schemes are also available to define the running coupling with the step scaling method. The gradient flow scheme is one of the applications of the gradient flow method, in which the gauge field is smeared with the so-called flow equation and the smeared gauge field has a nice perturbative property on the renormalizability [11][12][13]. In ref. [14], a renormalized coupling via the gradient flow in a finite size box with the periodic boundary condition has been introduced. However, Λ MS cannot be extracted from the coupling, since the coupling has a non-analytic expansion in α MS due to the zero-mode of the gauge field in the periodic boundary condition. To avoid the zero-mode problem, the twisted boundary condition has been introduced by Ramos [15]. The renormalized coupling defined in a finite box with the twisted boundary condition (the TGF scheme) has the normal one-loop relation to the MS scheme and is regularization independent. The running can be traced via the step scaling method on the lattice. The TGF running coupling for the pure SU(2) Yang-Mills theory has been evaluated using the step scaling method [15] and extended to two-color many flavor dynamical simulations [16]. The gradient flow coupling with the Schrödinger functional boundary condition is another scheme avoiding the zeromode problem and has been investigated in refs. [17][18][19] for the SU(3) gauge theories.
We extend Ramos's work [15] to the pure SU(3) Yang-Mills theory. In addition to this, we extract the Λ-parameter in the TGF scheme and convert it to the MS scheme. The ratio Λ MS /Λ TGF , which is usually evaluated using the perturbation theory, is not yet available at this time (but there is an ongoing study [20]). Since we already know Λ MS /Λ SF [10], actually what we have to estimate is the ratio Λ SF /Λ TGF . Therefore we estimate Λ SF /Λ TGF for the pure SU(3) gauge theory numerically with lattice simulations in this study. It should be noted that the analysis made in this paper is applicable to the gauge theories with dynamical fermions provided that the fermion representations and contents are compatible with the twisted boundary condition. This study is the first attempt to apply the TGF method for evaluating the Λ-parameter in the SU(3) gauge theories from the beginning to the end.
In this study we estimate Λ MS in terms of physical observables via the TGF method. Our strategy is summarized as follows: Here A phys is a physical observable with mass dimension and L max is an intermediate scale which connects the non-perturbative energy scale and the perturbative energy scale. In this paper, we consider the string tension √ σ and the Sommer scale 1/r 0 as the physical observable A phys . (Another reference scale can be considered, for example w 0 [21].) We will numerically calculate L max Λ TGF , L max A phys , and Λ SF /Λ TGF . L max Λ TGF is calculated with the step scaling method. In order to evaluate L max A phys , we employ data available from refs. [1,22] and ref. [23] for a √ σ and a/r 0 , respectively. We finally estimate Λ MS /A phys using eq. (1.3). We show that our estimates for Λ MS /A phys are compatible with the values previously obtained with other methods. This demonstrates the validity of our non-perturbative analysis. This paper is organized as follows. In the next section, we introduce the TGF method and explain how to calculate the TGF coupling briefly. Our strategy eq. (1.3) and the details of lattice simulations are explained in section 3. L max Λ TGF and L max A phys are presented in sections 4 and 5, respectively. Λ SF /Λ TGF and Λ MS /A phys are extracted in section 6. Finally we summarize this paper in the last section 7. Our preliminary result has been presented at the Lattice conference [24].

Twisted gradient flow coupling
We use the Wilson gauge action on a (L/a) 4 lattice with twisted boundary condition: Here U µ (n) is the SU(N C ) link variable with periodic boundary condition. We represent the twisted boundary condition by using the twist phase Z µν (n). In this work, we follow ref. [15] and put the twisted boundary condition in the x-y plane. The twist phase is defined as N C µ = 1, ν = 2, and n 1 = n 2 = 0, 1 otherwise, (2.2) in the case. The derivation of the action with the periodic variables (2.1) is given in appendix A. We first introduce link variables V µ (n, t) evolved with the gradient flow equation; where t, a fictitious time or so called flow time, is introduced. ∂ n,µ is the su(N C )-valued differential operator with respect to V µ (n, t).
The twisted gradient flow (TGF) coupling g 2 TGF (1/L) is defined as where E(t) is a energy density made of V µ (n, t). The explicit form of E(t) will be given later. The vacuum expectation value E(t) is a renormalized quantity at the scale 1/ √ 8t at any t > 0 [13]. In a finite volume system we can use the volume size L as the scale of the renormalization so we have set in eq. (2.4). The factor c is, in principle, a free parameter: a different choice of c gives a different renormalization scheme. Throughout this work we choose c = 0.3 for a reason we will state later. The normalization factor N −1 T (c, a/L) depends on the definition of the energy density on the lattice.
In this work, we employ the following definition for the energy density E(t); With this definition, the normalization factor N −1 T (c, a/L), which is defined so as to match g 2 TGF (1/L) with the bare coupling g 2 0 at the tree level of the perturbation theory, is The summation over P µ runs for µ = 1, 2 and for µ = 3, 4. The prime (′) symbol on the summation indicates the exclusion of the zero momentum contribution (P 1 , P 2 , P 3 , P 4 ) = (0, 0, 0, 0) from the sum. We employ c = 0.3 throughout this work. In general, a smaller value of c gives smaller statistical error. It causes, however, a larger lattice artifact. According to the previous works [15,17], c = 0.3 gives a good compromise between these two effects. This is the reason for our choice c = 0.3.

Overview of strategy and simulation details
Here we explain the strategy for evaluating eq.(1.3). We take the following steps.
1. We evaluate the discrete beta function B s (u) as a function of u = g 2 TGF (1/L). It is defined as where s is the scaling parameter. We extract this discrete beta function by taking the continuum limit of lattice discrete beta functions evaluated on several lattices. The details of the fitting and the analysis for the continuum limit will be explained in the next section.
2. We estimate L max Λ TGF using the discrete beta function evaluated in the previous step. By fixing the scale 1/L max implicitly through the value of the coupling u * = g 2 TGF (1/L max ), L max Λ TGF can be evaluated with Here we explicitly put c on the left-hand side, which is to use the same notation as eq. (2.5) for the scale setting. The TGF coupling at scale s n /L max is evaluated with the following recurrence equation (step scaling), For a sufficiently small value of u n = g 2 TGF (s n /L max ) we can safely use the two-loop approximation in eq. (3.2) to extract L max Λ TGF .
3. We relate the intermediate scale 1/L max to a hadronic scale A phys in the continuum limit. We employ two hadronic scales for the consistency check; the string tension √ σ and the Sommer scale r 0 . The lattice data of √ σ and r 0 are taken from refs. [1,22] and [23], respectively. To outline the procedure, let us assume that A phys has a mass dimension one for simplicity. We interpolate each of g 2 TGF (1/L, β) and aA phys (β) as a function of bare coupling β. By keeping the coupling constant g 2 TGF (1/L max , β * ) fixed to u * over several lattices L max /a * , we obtain the corresponding values of β * (here to show the connection between u * and the lattice spacing (or bare coupling), we use a * (or β * ) ). For each value of β * (thus a * /L max ) we have a pair of L max /a * and a * A phys (β * ). We then take the continuum limit of (L max /a * )(a * A phys ) as a function of a * /L max . 4. To convert Λ TGF to the Λ-parameter in the MS scheme, we need the ratio Λ MS /Λ TGF .
We split the ratio into two pieces: (Λ MS /Λ SF )(Λ SF /Λ TGF ). The value of the former factor is already known to be Λ MS /Λ SF = 0.48811(1) [10], but the latter is not known in the literature. We therefore calculate Λ SF /Λ TGF numerically via the one-loop relation between g 2 SF and g 2 TGF at the same renormalization scale 1/L. To obtain the one-loop relation, we calculate the couplings with lattice simulations in the weak coupling region.
5. Finally we combine the all pieces obtained above to have The TGF couplings on the lattice in the steps 1-3 explained above are evaluated on five lattices with L/a =12, 16, 18, 24 and 36. We use the heat-bath method introduced by Fabricius and Haan [25] to increase the acceptance ratio. We accumulate configurations as listed in table 7 in appendix B. Each configuration is separated by 100 sweeps. The TGF couplings, we computed, are listed in table 1, of which error is statistical one and estimated by taking the autocorrelation into account with the procedure proposed in ref. [26] 1 . We take several values for the bare coupling β = 6/g 2 0 on each lattice to take the continuum limit.
On the other hand, simulations in the weak coupling region have been done on four lattices L/a =8, 10, 12 and 16, with three values of the bare coupling β = 40, 60 and 80. We use the same plaquette gauge action with the O(a)-improvement boundary correction and the SF boundary condition [4,27] to calculate g 2 SF . The error of the coupling from these data is estimated with the Jackknife method after binning data into 10 bins. We execute O(10 6 )-O(10 7 ) sweeps for each parameter. The SF coupling is evaluated every sweep and the TGF coupling is evaluated every 100 sweeps. The error propagation of the statistical error on non-primary observables, such as the discrete beta function in the continuum limit, is estimated by a random re-sampling method. For the re-sampling, we assume the primary data in table 1 satisfies Gaussian distribution with the width of the measured statistical error.

TGF running coupling constant and L max Λ TGF
To extract the discrete beta function eq. (3.1), we take the continuum limit of the lattice discrete beta function defined by . (4.1) We use s = 3/2 as the scaling parameter. To take the continuum limit of eq. (4.1), the value of g 2 TGF (1/L, β) is kept fixed at g 2 TGF (1/L, β) = u as the renormalization condition irrespective of β. This implies that the physical length L is fixed. The lattice discrete beta function is evaluated using eq. (4.1) by substituting the data of table 1. We fit the lattice   discrete beta function with   5 Physical scale in terms of L max As described in section 3, the hadronic scales, the string tension √ σ and the Sommer scale r 0 , have to be determined in terms of L max . a √ σ and r 0 /a with the plaquette gauge action in large physical volumes have been determined at β ∈ [5.65, 6.515] in refs. [1,22] and β ∈ [5.70, 6.692] in ref. [23], respectively. To relate the intermediate scale L max /a and the physical scales aA phys (= a √ σ or a/r 0 ) at the same lattice cut-off "a", we need the bare coupling constant g 2 0 dependence (or β dependence) of g 2 TGF (1/L max , β) and aA phys (β). If the values of aA phys (β * ) at a fixed value g 2 TGF (a * /L max , β * ) = u * on several lattice sizes are obtained, we can take the continuum limit for L max A phys as follows:  Table 3. Fitted parameters for eq. (5.2) at each lattice size.
To take the continuum limit of the hadronic scale aA phys reliably, g 2 TGF (1/L max , β) should be precisely evaluated in the scaling region of aA phys on several lattice sizes L max /a with sufficiently large u * . This condition is satisfied with our data at L max /a = 12, 16 and 18, where the large enough TGF couplings g 2 TGF (1/L max , β) = u * and aA phys in the scaling region are available in the ranges β ∈ [6.11, 6.515] for a √ σ and β ∈ [6.11, 6.92] for a/r 0 , respectively. Therefore we can take any renormalization condition u * in this region and we employ several different values u * = 6.0, 6.1 . . . , 7.0 to see the consistency as stated in the previous section. Let us start with interpolation of g 2 TGF (1/L max , β), a √ σ and a/r 0 as functions of β separately in the following. Then we combine the interpolated results to take the continuum limit using eq. (5.1).
To interpolate g 2 TGF (1/L max , β), we fit the data at L max /a = 12, 16 and 18 in table 1 with the following interpolating function; We use all the data in β ∈ [6.1, 10.0] to stabilize the interpolation, while we need the interpolating formula only in the scaling region corresponding to u * we chose. Figure 3 shows the fit result, and table 3 shows the parameters obtained. It seems that χ 2 /DoF shown in table 3 are rather large, especially for L/a = 12. This is caused by using wider range than needed for the fitting. What we need is a smooth interpolating formula in the scaling region but not the fitting itself so we do not have to take the value of the χ 2 /DoF seriously. The scaling region of figure 3 is magnified in figure 4 showing a smooth interpolation of the fitting. Solving g 2 TGF (1/L max , β * ) = u * at each u * for β * using eq. (5.2), we obtain β * as shown in table 8 in appendix C.
Interpolating the data from refs. [1,22] for a √ σ as a function of β, we obtain with χ 2 /DoF ≃ 1.43. As plotted in figure 5, eq. (5.3) smoothly interpolates the data in the scaling region β ∈ [6.11, 6.515]. Substituting β * from table 8 into eq. (5.3), and multiplying L max /a * which corresponds to β * on it, we obtain L max √ σ at each u * . Table 9 in appendix C shows the values of L max √ σ before taking the continuum limit. The cut-off dependence of L max √ σ for each u * is shown in the left panel of figure 7. The values in the continuum limit are tabulated in the middle column of table 4.   in appendix C. The cut-off dependence and the values in the continuum limit are shown in the right panel of figure 7 and the right column in table 4, respectively.
6 Λ-parameter ratio Λ SF /Λ TGF and Λ MS To move from the TGF scheme to the MS scheme, we need the Λ-parameter ratio Λ MS /Λ TGF . Usually the ratio is calculated with the one-loop perturbation theory but the value is not yet available at the present time, while there is an ongoing project [20] of the perturbative   Table 4. L max √ σ and L max /r 0 for each u * in the continuum limit.
calculation. As we already know the ratio Λ MS /Λ SF [10], what we have to calculate is the ratio Λ SF /Λ TGF . Since both g 2 SF and g 2 TGF can be evaluated on the lattice with the same cut-off and with the renormalization scale (that is, a and L are the same), we can evaluate them with the Monte Carlo simulation on the lattice. We employ the two-loop formula [4,27], for the O(a)-improvement boundary correction in the SF simulations so that g 2 SF is O(a)improved at the two-loop level.
Let us denote the SF and TGF couplings at the gauge coupling β on a finite box (L/a) 4 by g 2 SF (a/L, β) and g 2 TGF (a/L, β), respectively. In a weak coupling region, these couplings  are related through We extract the value of c g (a/L) by investigating g 2 TGF (a/L, β) dependence of the ratio (6.2). Both couplings g 2 TGF and g 2 SF are numerically evaluated at β = 40, 60 and 80 on L/a = 8, 10, 12 and 16 lattices. Since the TGF scheme is automatically free from O(a) errors and g SF is O(a)-improved, the a/L dependence of c g (a/L) should be 3) The ratio of the Λ-parameters is defined by with c (0) g from the continuum limit of c g (a/L). In table 5 we list the TGF and SF couplings measured on each lattice size and each β. Figure 8 shows g 2 SF (a/L, β)/g 2 TGF (a/L, β) as a function of g 2 TGF (a/L, β). We fit the data linearly in g 2 TGF and the lines drawn in the figure are the fit results. Table 6 summarizes the fitted value of c g (a/L) for each L/a. In figure 9, we plot c g (a/L) as a function of (a/L) 2 . Fitting the data linearly in (a/L) 2 , we obtain  Table 6. The fit results for c g at each lattice.
We can now evaluate Λ MS according to our strategy eq. 0.48811(1) [10] and the results for L max Λ TGF , L max A phys , and Λ SF /Λ TGF (tables 2, 4 and eq. (6.6), respectively). Figures 10 and 11 show the renormalization condition u * dependence of Λ MS / √ σ and r 0 Λ MS , respectively. In these figures, the square symbols with error bar, which is statistical one, are our results and the dotted line is the average over our results with different u * . The dashed lines with gray band are from refs. [29] and [30] for comparison. We observe no renormalization condition dependence as expected. Our final estimates are

Summary
We have evaluated the Λ-parameter in the MS scheme for the pure SU(3) gauge theory via the twisted gradient flow method according to our strategy shown in (1.3). Our results are summarized in eqs. (6.7) and (6.8). To obtain the results we have determined the Λ-parameter ratio between the TGF scheme and the SF scheme with lattice simulations, which is a non-trivial step in our analysis. Having obtained sufficiently close values to the known ones in eqs. (6.7) and (6.8), we verified the ratio Λ SF /Λ TGF (6.6) determined with non-perturbative simulations. To further confirm the value of the ratio Λ SF /Λ TGF , it would be interesting to compare our ratio with the analytic one from the explicit perturbative calculation [20].

Acknowledgments
The numerical simulations have been done on the INSAM (Institute for Nonlinear Sciences and Applied Mathematics) cluster system at Hiroshima University. This work was partly supported by JSPS KAKENHI Grant Numbers 26400249 and 16K05326. I. K. is supported by MEXT as "Priority Issue on Post-K computer" (Elucidation of the Fundamental Laws and Evolution of the Universe) and JICFuS, by which K.-I. I. is also partially supported.

A Derivation of the action
In this appendix, we derive the action with periodic variables (2.1). We start from the following action in SU(N C ) defined on a (L/a) 4 ≡L 4 lattice with the twisted boundary condition on the x-y plane and periodic boundary condition in z and t directions: is a plaquette variable made of link variables V µ (n) with the twisted boundary condition: V µ (n +Lν) = V µ (n) (ν = 3, 4), (A.4) where N C × N C unitary matrix Γ ν (ν = 1, 2) is called twist matrix and satisfies Let us eliminate the link variables on n 1 = 0 or n 2 = 0 by using the variables on n i =L. The plaquette on n 1 = n 2 = 0 becomes Except for the overall factor ω, this is exactly the plaquette with periodic link variables U µ (n). Therefore we define link variables on n 1 = 0 and n 2 = 0 through the periodic boundary condition: U µ (0, n 2 , n 3 , n 4 ) ≡ U µ (L, n 2 , n 3 , n 4 ), (A.11) U µ (n 1 , 0, n 3 , n 4 ) ≡ U µ (n 1 ,L, n 3 , n 4 ), (A.12) U µ (0, 0, n 3 , n 4 ) ≡ U µ (L,L, n 3 , n 4 ). (A.13) Similar calculations show other plaquettes become those with U µ (n) without overall factor. Then, we finally obtain the action with periodic link variable U µ (n) Tr [Z µν (n)P µν [n; U ]] , (A.14) where Z µν (n) = Z * νµ (n) is given as Z µν (n) = ω * µ = 1, ν = 2, and n 1 = n 2 = 0, 1 otherwise. (A. 15) C Tables to evaluate L max /A phys In tables 8, 9 and 10 we collect values needed to evaluate L max /A phys in section 5.