High-precision quark masses and QCD coupling from n f =4 lattice QCD
Bipasha Chakraborty, C.T.H. Davies, G.C. Donald, R.J. Dowdall, B. Galloway, P. Knecht, J. Koponen, G.P. Lepage, C. McNeile
aa r X i v : . [ h e p - l a t ] D ec High-precision quark masses and QCD couplingfrom n f = 4 lattice QCD Bipasha Chakraborty, C. T. H. Davies, B. Galloway, P. Knecht, and J. Koponen
SUPA, School of Physics and Astronomy, University of Glasgow, Glasgow, G12 8QQ, UK
G. C. Donald
Institut f¨ur Theoretische Physik, Universit¨at Regensburg, 93040 Regensburg, Germany
R. J. Dowdall
DAMTP, University of Cambridge, Wilberforce Road, Cambridge, CB3 0WA, UK
G. P. Lepage ∗ Laboratory for Elementary-Particle Physics, Cornell University, Ithaca, NY 14853, USA
C. McNeile
School of Computing and Mathematics and Centre for Mathematical Science,Plymouth University, Plymouth PL4 8AA, United Kingdom (HPQCD Collaboration) (Dated: 4 December 2014)We present a new lattice QCD analysis of heavy-quark pseudoscalar-pseudoscalar correlators, using gluonconfigurations from the MILC collaboration that include vacuum polarization from u , d , s and c quarks ( n f =4 ). We extract new values for the QCD coupling and for the c quark’s MS mass: α MS ( M Z , n f = 5) =0 . and m c (3 GeV , n f = 4) = 0 . GeV. These agree well with our earlier simulationsusing n f = 3 sea quarks, vindicating the perturbative treatment of c quarks in that analysis. We also obtain anew nonperturbative result for the ratio of c and s quark masses: m c /m s = 11 . . This ratio implies m s (2 GeV , n f = 3) = 93 . MeV when it is combined with our new c mass. Combining m c /m s with ourearlier m b /m c gives m b /m s = 52 . , which is several standard deviations (but only 4%) away from theGeorgi-Jarlskop prediction from certain GUTs. Finally we obtain an n f = 4 estimate for m b /m c = 4 . which agrees well with our earlier n f = 3 result. The new ratio implies m b ( m b , n f = 5) = 4 . GeV.
PACS numbers: 11.15.Ha,12.38.Aw,12.38.Gc
I. INTRODUCTION
The precision of lattice QCD simulations has increaseddramatically over the past decade, with many calculationsnow delivering results with 1–2% errors or less. Such preci-sion requires increasingly accurate values for the fundamentalQCD parameters: the quark masses and the QCD coupling.Accurate QCD parameters are important for non-QCD phe-nomenology as well. For example, theoretical uncertaintiesin several of the most important Higgs branching fractionsare currently dominated by uncertainties in the heavy-quarkmasses (especially m b and m c ) and the QCD coupling [1].In this paper we present new lattice results for m c , m c /m s , m s , m b /m c , m b , and α s . In a previous paper [2] we ob-tained 0.6%-accurate results for the masses and couplingby comparing continuum perturbation theory with nonper-turbative lattice-QCD evaluations of current-current correla-tors for heavy-quark currents. Current-current correlatorsare particularly well suited to a perturbative analysis becausenon-perturbative effects are suppressed by four powers of Λ QCD / m h where m h is the heavy-quark mass. Our earlier ∗ [email protected] simulations treated u , d and s sea quarks nonperturbatively( n f = 3 ), while assuming that contributions from c and heav-ier quarks can be computed using perturbation theory. Herewe test the assumption that heavy-quark contributions are per-turbative by repeating our analysis with lattice simulationsthat treat the c quark nonperturbatively ( n f = 4 in the sim-ulation).In Section 2 we present our new n f = 4 lattice-QCD anal-ysis of current-current correlators, leading to new results forthe heavy-quark masses and the QCD coupling. We introducean improved procedure that gives smaller errors and simplifiesthe analysis. We also demonstrate how our Monte Carlo datacorrectly reproduce the running of the MS masses and cou-pling. In Section 3, we use the same simulations to calculatea new nonperturbative result for the ratio of the c to s quarkmasses, m c /m s . In Section 4, we use these simulations tocalculate the mass ratio m h /m c for heavy quarks with masses m h between m c and m b . We express the ratio as a functionof the heavy quark’s pseudoscalar mass m η h . We extrapo-late our result to m η h = m η b to obtain a new nonperturbativeestimate for m b /m c . In Section 5, we summarize our conclu-sions, derive new values for the s and b masses, and presentour thoughts about further work in this area. We also include,in Appendix A, a detailed discussion about how the couplingconstant, quark masses, and the lattice spacing depend uponsea-quark masses in our approach. Our current analysis in-cludes u/d sea-quark masses down to physical values, so weare able to analyze this in far more detail than before. Finally,Appendix B briefly summarizes n f = 4 results obtained usingour previous methods [2]. II. LATTICE RESULTS
Our new analysis follows our earlier work [2], but with asimpler and more accurate method for connecting current cor-relators to MS masses. In particular, this method allows usto determine the MS c mass at multiple scales, from correla-tors with different heavy-quark masses, providing a new testof our use of continuum perturbation theory. While the latticespacings are not as small as before, our new analysis treats c -quarks in the quark sea nonperturbatively. We also use thesubstantially more accurate HISQ discretization for the sea-quark action [3], in place of the ASQTAD discretization inour earlier analysis, and a more accurate method for settingthe lattice spacing. The gluon action is also improved overour earlier analysis, as it now includes O ( n f α s a ) correc-tions [4]. Our new results also have more statistics, and in-clude ensembles with u/d masses very close to the physicalvalue. A. Heavy-Quark Correlator Moments
As before, we compute (temporal) moments G n ≡ X t ( t/a ) n G ( t ) (1)of correlators formed from the pseudoscalar density operatorof a heavy quark, j ≡ ψ h γ ψ h : G ( t ) = a X x ( am h ) h | j ( x , t ) j (0 , | i . (2)Here m h is the heavy quark’s bare mass (from the latticeQCD lagrangian), a is the lattice spacing, time t is Euclideanand periodic with period T , and the sum over spatial posi-tions x sets the total three-momentum to zero. We again re-duce finite-lattice spacing, tuning and perturbative errors byreplacing the moments in our analysis with reduced moments: ˜ R n ≡ G /G (0)4 for n = 4 , m c (cid:16) G n /G (0) n (cid:17) / ( n − for n ≥ , (3)where G (0) n is the moment in lowest-order weak-coupling per-turbation theory using the lattice regulator, and m c is the baremass of the c quark.Low- n moments are dominated by short-distance physicsbecause the correlator is evaluated at zero total energy, whichis well below the threshold for on-shell hadronic states: thethreshold is at E threshold = m η h where2.9 GeV ≤ m η h < TABLE I. Perturbation theory coefficients for r n with n f = 4 seaquarks, where the heaviest sea quark has the same mass m h as thevalence quark (that is, the quark used to make the currents in thecurrent-current correlator). Coefficients are defined by r n = 1 + P j r nj α j MS ( µ ) where µ = m h ( µ ) . These coefficients are derivedin [6–10]. n r n r n r n . . − . . . − . . . . . . . for our range of masses m h . Furthermore, the moments areindependent of the ultraviolet cutoff when n ≥ . Apply-ing the Operator Product Expansion (OPE) to the product ofcurrents in the correlator, we can therefore write our n = 4 reduced moment in terms of continuum quantities, ˜ R → r ( α MS , µ ) d cond4 ( α MS , µ ) h α s G /π i eff (2 m h ) + ˜ d cond4 ( α MS , µ ) X q = u,d,s h m q ψ q ψ q i eff (2 m h ) + · · · , (5)in the continuum limit ( a → ). Here α MS is the MS couplingat scale µ , and m h is the MS h -quark mass. Heavy-quarkcondensates are absorbed into the gluon condensate [5]. Wewill retain terms only through the gluon condensate in whatfollows since its contribution is already very small and con-tributions from other condensates will be much smaller. Wediscuss the precise meaning of h α s G /π i eff below. Reducedmoments with n ≥ can be written: ˜ R n → r n ( α MS , µ ) m c ( µ ) (cid:26) d cond n ( α MS , µ ) h α s G /π i eff (2 m h ) + · · · (cid:27) , (6)where m c ( µ ) is the MS mass of the c quark. The contin-uum expressions for ˜ R n should agree with tuned lattice sim-ulations up to finite-lattice-spacing errors of O (( am h ) α s ) .The perturbative expansions for the coefficient functions r n are known through third order: see Table I and [6–10]. Theexpansions for d cond n are known through first order [11].Parameter µ sets the scale for m c and for α MS in r n . As inour previous paper, we take µ = 3 m h ( µ ) (7)in order to improve the convergence of perturbation theory. Infact, however, our method is almost completely independentof the choice of µ , by design. We can reexpress µ in terms ofthe MS mass of the c quark, µ = 3 m c ( µ ) m h m c , (8)since ratios of quark masses are regulator independent: thatis, m h m c = m h ( µ ) m c ( µ ) (9)up to a errors (for any µ ).Our reduced moments differ for n ≥ from our ear-lier work: here we multiply by /m c in Eq. (3) instead of m η h / m h . The ratio of G s in ˜ R n ≥ introduces a factorof m h /m h ( µ ) . This becomes /m c ( µ ) when multiplied by /m c (by Eq. (9)). Consequently we can use moments cal-culated with any heavy-quark mass m h to estimate the MS c mass (at µ = 3 m h ( µ ) ). Consistency among m c s comingfrom different m h values is an important test of the formal-ism.We could have used the bare mass of any quark, in place of m c , in Eq. (3). Then the n ≥ moments would give valuesfor the MS mass of that quark. Alternatively we could leavethe quark mass factor out, in which case these moments givethe factors Z m ( µ ) that convert any bare lattice quark mass intothe corresponding MS mass at scale µ . Heavy-quark current-current correlators, as used here, provide an alternative to RI-mom [12] and similar methods for determining both light andheavy quark masses.The new definition for the reduced moments simplifies ouranalysis since the variation of factor m c ( µ ) with µ is wellknown from perturbative QCD. The m η h dependence of theanalogous factor ( m η h / m h ) in the old analysis is unknown a priori , and so must be modeled in the fit. We analyzed ourdata using the old definitions; the results, which agree with theresults we find with the new methods, are described briefly inAppendix B. B. Lattice Simulations
To extract the coupling constant and c mass from simula-tions, we use the simulations to compute nonperturbative val-ues for the reduced moments ˜ R n with small n ≥ and a rangeof heavy-quark masses m h . We vary the lattice spacing, sowe can extrapolate to zero lattice spacing, and the sea-quarkmasses, so we can tune the masses to their physical values.The gluon-field ensembles we use come from the MILCcollaboration and include u , d , s , and c quarks in the quarksea [13, 14]. The parameters that characterize these ensemblesare given in Table II. The highly accurate HISQ discretiza-tion [3] is used here for both the sea quarks and the heavyquarks in the currents used to create the correlators. This dis-cretization was designed to minimize ( am h ) errors for large m h . Our previous work used HISQ quarks in the currents, buta less accurate discretization (ASQTAD) for the sea quarks.We also quote tuned values for the bare s and c quarkmasses in Table II. These are the quark masses that give the physical values for the η s and η c masses, as discussed in Ap-pendix A 1. This is the bare c mass we use in Eq. (3) for ˜ R n .In Table III we list our simulation results for the η h massand the reduced moments for various bare quark masses am h on various ensembles. Results from different values of am h on the same ensemble are correlated; we include these corre-lations in our analysis. The am η h values are computed fromBayesian fits of multi-state function X j =1 b j (cid:16) e − m j t + e − m j ( T − t ) (cid:17) (10)to the correlators G ( t ) for t ≥ , where T is the temporallength of the lattice [15]. The fitting errors are small for am η h and have minimal impact on our final results.The fractional errors in the ˜ R n for n ≥ are 20–40 timeslarger than those for ˜ R . This is because of the factor of /m tuned0 c used in Eq. (3) to define these moments. As men-tioned above, we could have used bare masses for otherquarks in this definition, to obtain values for their MS masses.Heavy-quark masses like m c , however, can usually be tunedmore accurately than light-quark masses, as discussed in Ap-pendix A. Masses for other quarks can be obtained from the c mass and nonperturbatively determined quark mass ratios,as we show for the s and b masses in the next two sections.As in our previous paper, we limit the maximum sizeof am h in our analysis: we require am h ≤ . . This keeps a errors smaller than 10%.We determine the lattice spacing by measuring the Wilsonflow parameter w /a on the lattice (Table II) [16]. From pre-vious simulations [17], we know that w = 0 . , (11)which we combine with our measured values of w /a to ob-tain the lattice spacing for each ensemble (Appendix A). Thisapproach is far more accurate than that used in our earlier pa-per, which relied upon the r parameter from the static-quarkpotential. C. Fitting Lattice Data
Our goal is to find values for α MS ( µ ) and m c ( µ ) that makethe theoretical results (from perturbation theory) for the re-duced moments ˜ R n (Eqs. (5–6)) agree with the nonperturba-tive results from our simulations. We do this by simultane-ously fitting results from all of our lattice spacings and quarkmasses for moments with ≤ n ≤ . To get good fits, wemust correct the continuum formulas in Eqs. (5–6) for sev-eral systematic errors in the simulation. We fit the lattice data TABLE II. Simulation parameters for the gluon ensembles used in this paper [13, 14], with lattice spacings of approximately 0.15, 0.12, 0.09and 0.06 fm, and various combinations of sea-quark masses. The parameters for each simulation are: the inverse lattice spacing in units of w = 0 . fm, the spatial L and temporal T lattice lengths, the number of gluon configurations N cf (each with multiple time sources),the bare sea-quark masses in lattice units ( am ℓ , am s , am c ), and the tuned bare s and c quark masses in GeV. The tuned s and c massesgives physical values for the η s and η c mesons, respectively. The ℓ mass is the average of the u and d masses, which are set equal in oursimulations. Z m ( µ ) is the ratio of the MS quark mass m q ( µ, n f = 4) to the corresponding bare (lattice) mass m q (see Section II D). Thelast two entries for each ensemble indicate the degree to which the sea-quark masses are detuned (see Appendix A).ensemble w /a L/a T /a N cf am ℓ am s am c m tuned0 s m tuned0 c Z m (3 GeV) δm sea uds /m s δm sea c /m c − .058(8)2 1.1272(7) 24 48 1000 0.00640 0.0640 0.828 0.0890(7) 1.130(4) 0.872(6) 0.046(14) − .050(8)3 1.1367(5) 36 48 1000 0.00235 0.0647 0.831 0.0885(7) 1.125(4) 0.876(5) − .048(13) − .034(8)4 1.3826(11) 24 64 300 0.01020 0.0509 0.635 0.0866(7) 1.057(3) 0.933(6) 0.236(16) − .044(8)5 1.4029(9) 32 64 300 0.00507 0.0507 0.628 0.0861(7) 1.051(3) 0.938(6) 0.067(14) − .035(8)6 1.4149(6) 48 64 200 0.00184 0.0507 0.628 0.0857(7) 1.047(3) 0.941(6) − .040(13) − .024(8)7 1.9330(20) 48 96 300 0.00363 0.0363 0.430 0.0823(9) 0.977(3) 1.009(6) 0.104(11) − .021(8)8 1.9518(7) 64 96 304 0.00120 0.0363 0.432 0.0818(7) 0.973(3) 1.013(6) − .011(13) − .003(8)9 2.8960(60) 48 144 333 0.00480 0.0240 0.286 0.0778(7) 0.912(3) 1.080(7) 0.365(19) 0.045(9)TABLE III. Simulations results for η h masses and reduced momentswith various bare heavy-quark masses am h and gluon ensembles(first column, see Table II). Only data for am h ≤ . are used infits to the correlators. am h am η h ˜ R ˜ R ˜ R ˜ R using the following corrected form: ˜ R n = ( for n = 41 /ξ m m c ( ξ α µ ) for n ≥ ) (12) × r n ( α MS ( ξ α µ ) , µ ) (13) × (cid:18) d cond n h α s G /π i (2 m h ) (cid:19) (14) × (cid:18) d h,cn m h − m c m h (cid:19) (15) + (cid:16) am η h . (cid:17) N X i =0 c i ( m η h , n ) (cid:16) am η h . (cid:17) i . (16)We use a Bayesian fit with priors for every fit parame-ter [15]. The priors are a priori estimates for the parametersbased upon theoretical expectations and previous experience,especially from our earlier, very similar n f = 3 analysis. Ineach case we test our choice of prior width against the Em-pirical Bayes criterion [15], which in effect uses fluctuationsin the data to suggest natural widths for priors. None of ourpriors is narrower than this optimal width, and most are wider,which leads to more conservative errors.We now explain each part of the lattice formula in turn.
1. Detuned Sea-quark Masses
The terms α MS ( ξ α µ ) and ξ m m h ( ξ α µ ) in ˜ R n are the MS coupling and heavy-quark mass for detuned sea-quark masses;see Eqs. (A9) and (A19) in Appendix A. Scale µ is chosen sothat µ = 3 ξ m m c ( ξ α µ ) m h m c = 3 m h ( µ, δm sea ) . (17)Scale factors ξ α and ξ m are defined in Appendix A, whichdiscusses how MS couplings and masses are affected by sea-quark masses. The coefficients g α , g m . . . in ξ α and ξ m aretreated as fit parameters, with priors taken from the output ofthe fits described in the appendix.The light sea-quark masses enter linearly in ξ α and ξ m , be-cause of (nonperturbative) chiral symmetry breaking. Quarkmass dependence also enters through the perturbation theoryfor the moments ( r n ), but is quadratic in the mass and there-fore negligible for light quarks. µ Dependence
The scale µ enters Eqs. (12)–(16) through the coupling con-stant α MS ( ξ α µ ) and the c mass m c ( ξ α µ ) . We parameterizethe coupling and mass in the fit by specifying their values at µ = 5 GeV with fit parameters α and m , α MS (5 GeV , n f = 4) = α m c (5 GeV , n f = 4) = m , (18)whose priors are α = 0 . ± . , m = 0 . ± . . (19)Our previous analysis gave . and . forthese parameters, so the priors are broad. The coupling andmass for other values of µ are obtained by integrating (nu-merically) their evolution equations from perturbative QCD,starting from the values at µ = 5 GeV: µ dα MS ( µ ) dµ = − β α ( µ ) − β α − β α − β α − β α , (20) d log m h ( µ ) d log µ = − γ α MS ( µ ) − γ α − γ α − γ α − γ α . (21)The first four coefficients on the right-hand-sides of theseequations are known from perturbation theory [18–21]. Ineach case, we treat the fifth coefficient as a fit parameter whoseprior’s width equals the root-mean-square average of the firstfour parameters: β = 0 ± σ β , γ = 0 ± σ γ . (22)Neither β nor γ has signficant impact on our final results.
3. Truncated Perturbation Theory
The Wilson coefficient function r n (Eq. (13)) has a pertur-bative expansion of the form r n ( α MS , µ ) ≡ N pth X j =1 r nj ( µ ) α j MS . (23)The perturbative coefficients r nj are known through third or-der, and are given for µ = m h ( µ ) in Table I. The lack of perturbative coefficients beyond third order isour largest single source of systematic error. Our data are suf-ficiently precise that higher-order terms are relevant. Further-more the relative importance of the higher-order terms varieswith quark mass, as α MS varies with µ = 3 m h . Thereforewe include the higher-order terms in our analysis with coeffi-cients that we fit to account for variations with quark mass. Asin our earlier analysis, we note that the known perturbative co-efficients are small and relatively uncorrelated from momentto moment and order to order for µ = m h , leading us to adoptfit priors r nj ( µ = m h ) = 0 ± (24)for the n > coefficients at µ = m h . We double the widthof these priors relative to our previous analysis because the fitsuggested that some higher-order coefficients are larger here(especially for n = 4 ).We set N pth = 15 terms in the expansion, although ourresults are essentially unchanged once 8 or more terms areincluded (or 5 with µ = m h ). As before we use renormal-ization group equations to express the coefficients r nj ( µ =3 m h ) in terms of the coefficients r nj ( µ = m h ) from Table Iand Eq. (24). This procedure generates (correlated) priorsfor the unknown coefficients at µ = 3 m h that account forrenormalization-group logarithms. The procedure makes ourresults largely independent of µ : our results change by lessthan a third of a standard deviation as µ is varied over theinterval m h ≤ µ ≤ m h .
4. Nonperturbative Effects; Finite-Volume Corrections
We use the Operator Product Expansion (OPE) in Eqs. (5–6) to separate short-distance from long-distance physics. Inprinciple, the perturbative coefficients in r n ( α MS , µ ) aboveshould have subtractions coming from the higher-order termsin the OPE expansion: r n → r n − d cond n h α s G /π i ( λ )pth (2 m h ) − · · · ! (25)where λ is a fixed cutoff scale in the perturbative regime, say λ = 1 GeV, and h α s G /π i ( λ )pth and d cond n are computed in per-turbation theory to the same order as r n . These subtractionscome from perturbative matching, and remove contributionsto r n due to low-momentum gluons ( q ≤ λ ), thereby also re-moving infrared renormalons order-by-order in perturbationtheory. The size of the subtraction depends upon the detaileddefinition of α s ( G ( λ ) ) . This procedure is completely unam-biguous given a specific definition for this operator, but wehave not included the subtraction in r n since it is negligiblefor any reasonable definition at our low orders of perturbationtheory. For example, a simple momentum-space cutoff, thatkeeps q < λ , gives [22] h α s G i ( λ )pth = 3 α s π λ , (26)which ranges from 0.001 to 0.019 GeV for λ s between500 Mev and 1 GeV. This would change r n by no more than0.1–0.4% at m h = m c and much less at our higher m h s.Not surprisingly, perturbative estimates of the condensatevalue (Eq. (26)) are similar in size to nonperturbative esti-mates of the condensate value. So it is simpler for us to com-bine the subtraction in Eq. (25) with the condensate itself toform an effective condensate value [23]: h α s G i eff ≡ h α s G i ( λ ) − h α s G i ( λ )pth (27)In our fits we take h α s G i eff as a fit parameter with prior h α s G i eff = 0 . ± . , (28)and we approximate m h ≈ m η h / . in the condensate cor-rection (because m b ( m b ) ≈ m η b / . ). Our results are com-pletely unchanged if the width of this prior is ten times larger.In either case we obtain a value for the effective condensate oforder . with errors of a similar size. This is completelyconsistent with expectations, and it reduces condensate con-tributions to the moments to 0.01–0.05% at m h = m c , andmuch less at higher m h — negligible at our level of precision.This procedure is sensible at our level of precision. As pre-cision increases, however, there is a point where it becomesimportant to remove renormalon corrections from the coeffi-cients in r n . Otherwise j ! factors in j th order, coming frominfrared renormalons, cause perturbation theory to diverge. Asimple analysis [24] indicates that perturbation theory startsto diverge at order j ∼ / ( β α MS ) , which is around th or-der for our analysis. Consequently we expect the impact ofinfrared renormalons to be negligible at rd order.Perturbation theory is not the whole story even if in-frared renormalons are removed. The OPE separates short-distances from long-distances, but the short-distance coef-ficients r n , d cond n . . . have nonperturbative contributions, forexample, from small instantons [22]. It is also possiblethat the OPE is an asymptotic expansion and does not con-verge ultimately, although recent results suggest it might con-verge [25, 26]. Whatever the case, such effects are expectedto appear at even higher orders than infrared renormalons, andso are completely negligible at our level of precision.Condensates, renormalons, small instantons, etc. afflict allperturbative analyses at some level of precision. Our analysisis particularly insensitive to such effects because the leadingnonperturbative contributions are suppressed by four powersof Λ QCD / (2 m h ) .Note finally that the coefficient functions, being short-distance, are insensitive to errors caused by the finite volumeof the lattice. While the finite volume can affect the valueof h α s G i eff , the impact on our results is negligible since thecondensate itself is negligible. We verified this by recalcu-lating the reduced moments for emsemble 5 in Table II withspatial lattice sizes of L/a = 24 and 40 (ensemble 5 uses 32).The moments for different volumes agree to within statisticalerrors of order 0.01%. The same is true for the measured val-ues of m η c from these ensembles; finite volume effects willbe smaller still for m η h . m h − m c Correction
Our results are also affected by the difference between the c mass m c used in the sea, and the mass of the heavy quark m h used to make the currents in the current-current correla-tor. The perturbative calculations we use assume m c = m h ,but we want to study a range of m h values with fixed m c .The correction enters in O ( α s ) , is quadratic in the mass dif-ference for small differences, and goes to a (small) constantas m h → ∞ . Therefore we correct for it using (Eq. (15)) ˜ R n → ˜ R n (cid:18) d h,cn m h − m c m h (cid:19) (29)where h n is a fit parameter with a prior of ± . . Thewidth 0.03 is ten times larger than the correct value (from per-turbation theory) in the m h → ∞ limit. It is twice as wideas the width indicated by the Empirical Bayes criterion [15].We also tried fits where d h,cn was replaced by a spline func-tion of m η h . These give similar results but with larger errors(especially for α MS ).
6. Finite Lattice Spacing Errors
The final modification in our formula for ˜ R n corrects forerrors caused by the finite lattice spacings used in the simula-tions. We write ˜ R n → ˜ R n + δ ˜ R n (30)where δ ˜ R n ≡ (cid:16) am η h . (cid:17) N X i =0 c ( n ) i ( m η h ) (cid:16) am η h . (cid:17) i (31)and again m η h / . is a proxy for the quark mass. We pa-rameterize the m η h dependence of the c ( n ) i ( m η h ) using cubicsplines with knots, at m knots ≡ { } GeV , (32)that come from the analysis in Section IV. We set c ( n ) i ( m ) = c ( n )0 i + δc ( n ) i ( m ) (33)with the following fit parameters and priors: c ( n )0 i = 0 ± /nδc ( n ) i ( m ) = 0 ± . /n m ∈ m knots δc ( n ) ′ i ( m ) = 0 ± . /n m = 2 . . (34)These priors are again conservative since the Empirical Bayescriterion [15] suggests priors that are half as wide. We take N = 20 but our results are insensitive to any N ≥ . D. n f = 4 Lattice Results
We fit all of the reduced moments from our simulationdata — with lattice spacings from 0.12 fm to 0.06 fm, and n = 4 , 6, 8 and 10 in Table III — simultaneously to for-mula (12–16) by adjusting fit parameters described in the pre-vious sections. The fit is excellent with a χ per degree offreedom of 0.51 for 92 pieces of data ( p -value is 1.0).The fit has two key physics outputs. One is a new result forthe running coupling constant: α MS (5 GeV , n f = 4) = 0 . . (35)To compare with our old determination and other determi-nations, we use perturbation theory to add b quarks to thesea [27], with m b ( m b ) = 4 . GeV [2], and evolve tothe Z mass (91.19 GeV) to get α MS ( M Z , n f = 5) = 0 . . (36)This agrees well with . from our n f = 3 analysis [2].It also agrees well with the current world average 0.1185(6)from the Particle Data Group [28].The second important physics output is the c quark’s mass,whose value at µ = 5 GeV is a fit parameter: m c ( µ, n f = 4) = . µ = 5 GeV0 . µ = 3 GeV1 . µ = m c ( µ ) , (37)where we have used Eq. (21) to evolve our result to otherscales for comparison with other determinations. Theseagain agree well with our previous n f = 3 analysis [2],which gave 0.986(6) GeV for the mass at 3 GeV. The errorsfor m c (3 GeV) and α MS ( M Z ) are correlated, with correla-tion coefficient 0.19.We use our result from m c to calculate the mass renormal-ization factors Z m ( µ ) ≡ m c ( µ ) m c (38)that relate MS masses to bare lattice masses for each config-uration. These factors can be used to convert the bare massfor any quark to its MS equivalent. We tabulate these results,with µ = 3 GeV, for our configurations in Table II. These Z m values are much more accurate than can be obtained fromorder α s lattice QCD perturbation theory [29], but they agreequalitatively and suggest that higher-order corrections fromlattice perturbation theory are small.Our results confirm that a perturbative treatment of c quarksin the sea, as in our previous paper, is correct, at least to ourcurrent level of precision.Our result at µ = m c has a larger error because α MS inthe mass evolution equation (Eq. (21)) becomes fairly largeat that scale ( α MS ≈ . ) and quite sensitive to uncertaintiesin its value. We use the coupling from our fit for this evolu-tion. Were we instead to use the Particle Data Group’s (moreaccurate) α MS , our value for m c ( m c ) would be m c ( m c , n f = 4) = 1 . . (39) . . . . m c ( m h ) n = . . . . m c ( m h ) n = . . . . . m h / m c . . . . m c ( m h ) n = FIG. 1. The c quark mass m c ( µ = 3 m h ) as determined from mo-ments with heavy-quark masses ranging from m c to . m c . Thedata points show results obtained by substituting nonperturbativesimulation values for ˜ R n into Eq. (40), after correcting for mistun-ings of the sea-quark masses (using the fit). Errors are about thesize of the plot symbols, or smaller. Results are shown for threelattices spacings: 0.12 fm (green points, through m h /m c = 1 . ),0.09 fm (blue points, through m h /m c = 1 . ), and 0.06 fm (redpoints, through m h /m c = 2 . ). The dotted lines show our fits tothese data points. The gray band shows the values expected from ourbest-value m c (5 GeV) = 0 . GeV evolved perturbatively tothe other scales.
In any case, it is probably better to avoid such low scales, ifpossible.Note that our c mass comes from moments whose heavy-quark mass varies from m h = m c to m h = 3 m c . Each (non-perturbative) ˜ R n with n ≥ , for each heavy-quark mass m h ,gives an independent estimate of the c mass: m c (3 m h ) = r n ( α MS (3 m h ) , µ = 3 m h )˜ R n . (40)The extent to which these estimates agree with each other isshown in Figure 1, where the nonperturbative results (datapoints) are compared with our best-fit result for m c (5 GeV) evolved perturbatively to other scales using Eq. (21) (grayband). As expected, finite a errors are larger for smaller val-ues of n and larger values of m h [2, 30]. Taking account ofthese errors, agreement between different determinations ofthe mass is excellent.The dominant sources of error for our results are listed inTable IV. The most important systematics are due to the trun-cation of perturbation theory and our extrapolation to a = 0 .As in our previous analysis, the a extrapolations are not TABLE IV. Error budget [31] for the c mass, QCD coupling, andthe ratios of quark masses m c /m s and m b /m c from the n f = 4 simulations described in this paper. Each uncertainty is given as apercentage of the final value. The different uncertainties are added inquadrature to give the total uncertainty. Only sources of uncertaintylarger than 0.05% have been listed. m c (3) α MS ( M Z ) m c /m s m b /m c Perturbation theory 0.3 0.5 0.0 0.0Statistical errors 0.2 0.2 0.3 0.3 a → δm sea uds → δm sea c → m h = m c (Eq. (15)) 0.1 0.1 0.0 0.0Uncertainty in w , w /a α prior 0.0 0.1 0.0 0.0Uncertainty in m η s m h /m c → m b /m c δm η c : electromag., annih. 0.1 0.0 0.1 0.1 δm η b : electromag., annih. 0.0 0.0 0.0 0.1Total: 0.64% 0.63% 0.55% 1.20% . . . . . a (GeV − )0 . . . . . . ˜ R n ˜ R ˜ R ˜ R ˜ R FIG. 2. Lattice-spacing dependence of reduced moments ˜ R n for η h masses within 5% of m η c , and n = 4 , 6, 8, 10. The dashedlines show our fit, and the points at a = 0 are the continuum extrap-olations of the lattice data. large, as is clear from Figure 1 and also Figure 2. Also the de-pendence of our results on the light sea-quark masses is quitesmall and independent of the lattice spacing, as illustrated byFigure 3.Our results change by σ/ if we fit only the n = 4 and 6moments, but the errors are 35% larger. Leaving out n = 4 ,instead, leaves the c mass almost unchanged, but increases theerror in the coupling by 60% (with the same central value).We limit our analysis to heavy quark masses with am h ≤ . , as in our previous analysis. Reducing that limit to . , forexample, has no impact on the central values of results andincreases our errors only slightly (less than 10%).We tested the reliability of our error estimates for the per-turbation theory by refitting our data using only a subset ofthe known perturbative coefficients. The results are presentedin Fig. 4, which shows values for m c (3 GeV) and α MS ( M Z ) . . . . . ˜ R . . . . . ˜ R . . . . . ˜ R . . . . δ m sea uds / m s . . . . . ˜ R FIG. 3. Light sea-quark mass dependence of reduced moments ˜ R n for m h = m c , and n = 4 , 6, 8, 10. Results are shown for our twocoarsest lattices: a = 0 . fm (three points in blue) and a = 0 . fm(two points in red). The dashed lines show the corresponding resultsfrom our fit. Note that the slopes of the lines are independent of thelattice spacing, as expected. from fits that treat perturbative coefficients beyond order N as fit parameters, with priors as in Eq. (24). Results from dif-ferent orders agree with each other, providing evidence thatour estimates of truncation errors are reliable. This plot alsoshows the steady convergence of perturbation theory as addi-tional orders are added.As a further test of perturbation theory, we refit our nonper-turbative data treating the leading perturbative coefficients, γ and β , in the evolution equations for the mass (Eq. (21)) andcoupling (Eq. (20)) as fit parameters with priors of ± . Thefit gives γ = 0 . β = 0 . , (41)in good agreement with the exact results of . and . ,respectively. So our nonperturbative results for the correlatorsshow clear evidence for the evolution of m c ( µ ) and α MS ( µ ) as µ = 3 m h varies from m c to m c . . . . . . m c ( G e V ) N . . . . . α M S ( M Z ) FIG. 4. Results for the MS c mass and coupling from n f = 4 fitsthat treat perturbative coefficients beyond order N as fit parameters,with priors specified by Eq. (24). The gray bands and dashed linesindicate the means and standard deviations of our final results, whichcorrespond to N = 3 . . . . . . . . . . ( am c ) . . . . . . . m c / m s FIG. 5. The ratio of the c and s quark masses as a function of thesquared lattice spacing (in units of the bare c mass). The data comefrom simulations at lattice spacings of 0.15, 0.12, 0.09 and 0.06 fm,after tuning the s and c masses to reproduce physical values for the η s and η c masses on each ensemble. The errors for the data points arehighly correlated, as they come primarily from uncertainties in w , m η s , and m η c . The red dashed line shows our fit, which has a χ perdegree of freedom of 0.21 for 9 degrees of freedom ( p -value of 0.99).The black dashed line and gray band show the mean value and stan-dard deviation for our result extrapolated to zero lattice spacing. III. m c /m s FROM n f = 4 As discussed above (Section II A), we can use lattice QCDto extract ratios of MS quark masses completely nonperturba-tively [32], since ratios of quark masses are scheme and scale independent: for example, m c m s (cid:12)(cid:12)(cid:12)(cid:12) lat = m c ( µ, n f ) m s ( µ, n f ) (cid:12)(cid:12)(cid:12)(cid:12) MS + O (( am c ) α s ) . (42)While ratios of light-quark masses can be obtained from chiralperturbation theory, only lattice QCD can produce nonpertur-bative ratios involving heavy quarks. These ratios are veryuseful for checking mass determinations that rely upon per-turbation theory, as illustrated in [2]. They also allow us toleverage precise values of light-quark masses from very accu-rately determined heavy-quark masses.In [32] we used nonperturbative simulations, with n f = 3 sea quarks, to determine the s quark’s mass from the c quark’smass and the ratio m c /m s . We repeat that analysis here, butnow for n f = 4 sea quarks, using the tuned values of the bare s and c masses for each of our lattice ensembles: am tuned0 s and am tuned0 c in Table II, respectively. We expect am tuned0 c am tuned0 s = m c m s h m δm sea uds m s + h a ,m δm sea uds m s (cid:18) m c π/a (cid:19) + h α s ( π/a ) (cid:18) m c π/a (cid:19) + N a X j =2 h j (cid:18) m c π/a (cid:19) j , (43)where again we ignore δm sea c and δm dependence since theyare negligible. We fit the data from Table II using this formulawith the following fit parameters and priors: h m = 0 ± . , h a ,m = 0 ± . , (44) h = 0 ± , h j = 0 ± j > . (45)The extrapolated value m c /m s is also a fit parameter. We set N a = 5 , but get identical results for any N a ≥ .The result of this fit is presented in Fig. 5, which showsthe a dependence of the lattice results. The sensitivity of ournew results to a is about half what we saw in our previousanalysis. Our new fit is excellent and gives a final result forthe mass ratio of: m c ( µ, n f ) m s ( µ, n f ) = 11 . . (46)The leading sources of error in this result are listed in Ta-ble IV. These are dominated by statistical errors and uncer-tainty in the η s mass. Many other potential sources of error,such as uncertainties in the lattice spacing, largely cancel inthe ratio.Note that the discussion in Appendix A and Eq. (A19),in particular, imply that the leading effect of mistuned sea-quark masses cancels in ratios of quark masses. This is sub-stantiated by our fit which makes parameter h m negligiblysmall ( − . ). Setting h m = 0 shifts our result for m c /m s by only σ/ .Our result is a little more than a standard deviation lowerthan the recent result, . (cid:0) +59 − (cid:1) , computed by the Fer-milab/MILC collaboration (using many of the same configu-rations we use) [33]. Our analysis uses a different scheme for0 m h / m c n f = m η h (GeV)1234 m h / m c n f = FIG. 6. The ratio of the h and c quark masses as a function of themass of hh pseudoscalar meson mass. The data come from simu-lations at lattice spacings of 0.15, 0.12, 0.09 and 0.06 fm; the datapoints are colored magenta, blue, green, and red, respectively. Thegray band and dashed line in the top panel show function Eq. (47)with the best fit parameters, extrapolated to zero lattice spacingand the correct sea-quark masses. The bottom panel compares the n f = 4 data with extrapolated results obtained in [2] from current-current correlators in n f = 3 simulations. tuning the lattice spacing and quark masses, which leads tothe lack of sea-quark mass dependence in m c /m s discussedjust above. The absence of sea-mass dependence is apparentfrom Fig. 5, where the clusters of data points correspond to en-sembles with the same bare lattice coupling but different sea-quark masses. This figure can be compared with Fig. 6 in [33],which shows much larger sea-mass dependence. Both ap-proaches should agree when extrapolated to zero lattice spac-ing and the physical sea-quark masses. IV. m h /m c FROM m η h An analysis similar to that in the previous section allows usto relate heavy-quark masses m h to the hh pseudoscalar mass m η h with data from Table III. This can be used, for example,to estimate the b mass by extrapolating to m η b .Here we fit the lattice mass ratios m h /m tuned0 c to the fol-lowing function of m η h from the simulation: m h m c = m η h m η c N X n =0 f n ( m η h ) (cid:16) am η h (cid:17) n + f sea ( η h ) m η h m η c δm sea uds m s (cid:16) am η h (cid:17) (47)where N = 20 , although any N > gives the same result.Here f n ( m η h ) and f sea ( m η h ) are cubic splines with knots at m knots = { . , . , . , . } GeV . (48)The maximum and minimum knots correspond to the maxi-mum and minimum values of m η h , while the locations of the internal knots were obtained by treating those locations as fitparameters. Each f is parameterized by f ( m ) = f + δf ( m ) (49)and fit parameters f = 0 ± δf ( m ) = 0 ± . m ∈ m knots δf ′ ( m ) = 0 . ± . m = 2 . . (50)We reduce the priors for the leading a errors by a factorof / since these errors are suppressed by α s in the HISQ dis-cretization. The choice of priors for the spline parameters ismotivated by results from [2] (see Figure 4 in that paper).The fit is excellent with a χ per degree of freedom of 0.44for 29 pieces of data: see the top panel in Figure 6. Finitelattice spacing errors are much smaller for this quantity thanfor the moments, and it is again largely independent of mis-tunings in the sea-quark masses. Extrapolating to m η b gives m b /m c = 4 . (51)which agrees with our n f = 3 result of . , but withlarger errors [2]. Our new n f = 4 data go down to lat-tice spacings of . fm; our earlier analysis also had resultsat . fm.The bottom panel of Figure 6 compares our new n f =4 data with n f = 3 results obtained from fits to the current-current correlators [2]. The agreement is excellent, showingagain that n f = 3 and n f = 4 are consistent with each other. V. CONCLUSIONS AND OUTLOOK
The initial extractions of quark masses from heavy-quarkcurrent-current correlators relied upon experimental data from ee annihilation [34, 35]. Our analysis here, like the two thatpreceded it [2, 30], replaces experimental data with nonper-turbative results from tuned lattice simulations.Lattice simulations offer several advantages over experi-ment for this kind of calculation [1]. For one thing, simu-lations are easier to instrument than experiments and muchmore flexible. Thus we can generate lattice “data” not justfor vector-current correlators, but for any heavy-quark cur-rent or density; we optimize our simulations by using thepseudoscalar density instead of the vector current. Experi-ment provides results for only two heavy-quark masses — m c and m b — but we can produce lattice data for a whole rangeof masses between m c and m b . This means that α MS ( µ ) varies continuously, by almost a factor of two, in our analysissince µ ∝ m h . Here we use this variation to estimate andbound uncalculated terms in perturbation theory, providingmuch more reliable estimates of perturbative errors than thestandard procedure of replacing µ by µ/ and µ . (Our anal-ysis is essentially independent of µ .) Nonperturbative contri-butions are also strongly dependent upon m h , and thereforemore readily bound if a range of masses is available; they arenegligible in our analysis.1 . . . . . m c / m s HPQCD 0910.3102ETMC 1010.3659ETMC 1403.4504MILC 1407.3772HPQCD this paperDurr 1108.1650u,d,s,c seau,d,s seau,d sea
FIG. 7. Lattice QCD determinations of the ratio of the c and s quarks’masses. The ratios come from this paper and references [32, 33, 36–38]. The gray band is the weighted average of the three n f = 4 results: . . In this paper, we have redone our earlier n f = 3 analysis [2]using simulations with n f = 4 sea quarks: u , d , s and c . Ournew results, m c (3 GeV , n f = 4) = 0 . (52) α MS ( M Z , n f = 5) = 0 . , (53)agree well with our earlier results of . GeV and . , suggesting that contributions from c quarks inthe sea are reliably estimated using perturbation theory (asexpected). Our c mass is about . σ lower than the re-cent result from the ETMC collaboration, also using n f =4 simulations but with a different method [36]: they get m c ( m c ) = 1 . GeV, compared with our n f = 4 re-sult of 1.2715(95) GeV.We updated our earlier n f = 3 analysis [32] of the ra-tio m c /m s of quark masses using our n f = 4 data. Thisis a relatively simple analysis of data from Table II. Our newvalue is: m c ( µ, n f ) m s ( µ, n f ) = 11 . . (54)It agrees well with our previous result . , but is muchmore accurate. We compare our new result with others inFig. 7.We obtain a new estimate for the s mass by combining ournew result for m c /m s with our new estimate of the c mass(Eq. (52), converted from n f = 4 ): m s ( µ, n f = 3) = ( . µ = 2 GeV84 . µ = 3 GeV . (55)Values for m s ( µ, n f = 4) are smaller by about 0.2 MeV. Ournew result agrees with our previous analysis and also with other recent n f = 3 or 4 analyses: m s (2 GeV) = . .
5) MeV
HPQCD [32] , . .
1) MeV
ETMC [36], . .
9) MeV
Durr et al [39], m s (3 GeV) = 83 . .
0) MeV
RBC/UKQCD [40] . (56)Finally, we have also updated our previous ( n f = 3 ) non-perturbative analysis of m b /m c using our new n f = 4 data.We obtain: m b ( µ, n f ) m c ( µ, n f ) = 4 . , (57)which agrees with our previous result of 4.51(4) [2]. Combin-ing this result with our new value for m c (Eq. (52)) gives m b ( m b , n f = 5) = 4 . . (58)This again agrees with our earlier result of 4.164(23) GeV, butwith larger errors. We can also multiply our results for m b /m c and m c /m s to obtain m b ( µ, n f ) m s ( µ, n f ) = 52 . . (59)This is almost four standard deviations (but only 4%) awayfrom the result predicted by the Georgi-Jarlskog relation-ship [41] for certain classes of grand unified theory: theGeorgi-Jarlskog relationship says that m b /m s should equal m τ /m µ = 50 . .The prospects for improving our results over the nextdecade are good. Detailed meta-simulations, described in [1],indicate that errors from our analysis can be pushed below0.25% by a combination of higher-order perturbation the-ory, and, especially, smaller lattice spacings (0.045, 0.03and 0.023 fm) — both improvements that are quite feasibleover a decade [1]. There are also many other promising ap-proaches within lattice QCD. Several exist already for extract-ing the QCD coupling: see, for example, [42–47]. One canalso use simulations of other renormalized quantities, such asthe m h ψ h γ ψ vertex function, to compute quark masses [12].Small lattice spacings are particularly important for the b mass, because lattice spacing errors are typically of or-der ( am b ) . One approach is to use highly-improved relativis-tic actions for the b quarks, like the HISQ action used here. Asshown in [3], all but one of the O ( a, a ) operators that arise inthe Symanzik improvement of a quark action are suppressedby extra factors of the heavy-quark velocity: factors of ( v/c ) for mesons made of heavy quarks, and v/c for mesons madeof a combination of heavy and light quarks. The one opera-tor that does not have extra suppression is P µ ψγ µ ( D µ ) ψ ,which violates Lorentz invariance and so is easily tuned non-perturbatively using the meson dispersion relation. This is thestrategy adopted in the HISQ discretization we use here. Theextra factors of v/c suppress ( am b ) errors by an extra or-der of magnitude, beyond the suppression, by a power of α s ,coming from tree-level corrections for a errors in HISQ.2 ( am b ) errors can be avoided completely by using effec-tive field theories like NRQCD [48] or the Fermilab formal-ism [49] for b dynamics. Such approaches should be suf-ficiently accurate provided they are corrected to sufficientlyhigh order in ( v b /c ) . Our recent NRQCD analysis of m b ,using current-current correlators, is encouraging [50].Overall the prospects are excellent for continued improve-ment. ACKNOWLEDGMENTS
We are grateful to the MILC collaboration for the use oftheir gauge configurations and code. We thank S. King andD. Toussaint for useful conversations. Our calculations weredone on the Darwin Supercomputer as part of STFC’s DiRACfacility jointly funded by STFC, BIS and the Universities ofCambridge and Glasgow. This work was funded by STFC,the Royal Society, the Wolfson Foundation and the NationalScience Foundation.
Appendix A: Sea-Quark Mass Dependence
In this appendix we discuss the dependence of the MS cou-pling and heavy-quark masses on the sea-quark masses. Wevary the u/d sea-quark mass in our simulations to help us as-sess systematic errors associated with tuning that mass. In ad-dition, the precision with which the s and c sea-quark masseshave been tuned varies by several percent over the various en-sembles we use. These detunings shift the MS coupling andmasses. We need to understand how they are shifted in or-der to extract results for α MS and m h with physical sea-quarkmasses.It is essential when discussing detuned sea-quark masses tobe specific about what is held fixed as the quark masses areshifted from their physical values. An obvious choice is tofix both the lattice spacing a and the bare coupling α lat in thelattice lagrangian, while varying the quark masses. We findit more convenient, however, to explore a slightly differentmanifold in theory space by fixing α lat and the value of theWilson-flow parameter w .Lattice simulations are done for particular values of the barecoupling constant (and bare quark masses), but with all di-mensional quantities expressed in units of the lattice spacing( lattice units ). This removes explicit dependence on the lat-tice spacing from the simulation, so we can run the simulationwithout knowing the lattice spacing. To extract physics, how-ever, we must determine the lattice spacing (from the sim-ulation) and convert all simulation results from lattice unitsto physical units. In our simulations, we calculate the latticespacing by measuring the value of a/w in the simulation, andmultiplying it by the known value of w for physical sea-quarkmasses (that is, 0.1715(9) fm). As a result the lattice spacingbecomes (weakly) dependent upon the sea-quark masses since w is affected by sea quarks.This procedure is convenient because the lattice spacing fora given ensemble is determined using information from only that ensemble, thereby decoupling the analyses of differentensembles to a considerable extent. As we discuss below thereis an added benefit when vacuum polarization from c (or heav-ier) quarks is included in the simulation, as we do here: heavyquarks automatically decouple from low-energy physics (like w [51]). With our procedure, physical quantities that probeenergy scales smaller than m c — that is, almost everythingstudied with lattice QCD today — are essentially independentof m c , which means that they are completely unaffected bytuning errors in m c . This would not be the case if we fixed thelattice spacing instead of w , since it is small variations in thelattice spacing that correct for mistuning in m c .It is also very convenient that we set the lattice spacing us-ing a flavor singlet quantity. Because w is a flavor singlet, theleading sea-mass dependence induced in the lattice spacing isanalytic (linear) in the quark mass and small; in particular,there are no chiral logarithms [52]. One consequence is thatleading-order chiral perturbation theory for physical quanti-ties ( f π , f D s . . . ) is unchanged from standard treatments ex-cept for shifts (that are easily accommodated) in the coeffi-cients of certain analytic terms.In this appendix we show how the MS coupling and heavy-quark mass depend upon the sea-quark masses in our simu-lations. This dependence implies sea-quark mass dependencein the lattice spacing and the heavy quark’s bare mass, whichwe then use to determine some of the parameters involved.Finally we review heavy-quark decoupling, and estimate theparameters for c -mass dependence using first-order perturba-tion theory.
1. Tuning Bare Quark Masses
We define tuned values for the bare c and s masses on eachensemble by adjusting those masses to give physical values insimulations for the η c and η s masses. The tuned values arelisted in Table II.The current experimental value for the η c massis 2.9836(7) GeV [28]. In our analysis, we remove electro-magnetic corrections from this value, and adjust its errorto account for cc annihilation, since neither effect is in oursimulations [53, 54]. We use: m phys η c = 2 . . (A1)We compute the tuned c mass m tuned0 c by linear interpolationusing η h masses from the simulation (Table III) for heavy-quark masses m h in the vicinity of m c . In a few cases wehave results for only a single value of m h ; then we computethe tuned c mass using estimates of dm η c /dm c from otherensembles with (almost) the same lattice spacing.Note that the uncertainty in m tuned0 c is usually smaller thanthat in am tuned0 c . This is a peculiar feature of heavy-quarkmasses in lattice simulations (see, for example, [55]). It fol-lows from the formula for the linear interpolation that definesthe tuned mass in terms of a nearby mass: m tuned0 c = ( am c ) a − + dm c dm η c (cid:0) m phys η c − ( am η c ) a − (cid:1) (A2)3where am η c is the simulation result for the η c mass (in latticeunits) when the c quark has mass am c . Here dm c /dm η c isobtained from simulation results for a few nearby c masses.The uncertainty in a − is usually larger than the uncertaintiesin the other lattice quantities, but here a − is multiplied by ( am c ) − ( am η c ) dm c dm η c (A3)which would vanish if m η c = 2 m c . This cancellationis only partial for real masses, but it doesn’t occur at allif Eq. (A2) is multiplied on both sides by a to give a for-mula for am tuned0 c . As a result, fractional errors are roughly × smaller for m tuned0 c .The η s is an ss pseudoscalar particle where the valencequarks are (artificially) not allowed to annihilate; its physi-cal mass is determined in lattice simulations from the massesof the pion and kaon [17]: m phys η s = 0 . (A4)This mass is defined for use in lattice simulations and needs nofurther corrections for electromagnetism. We tune the s massby simulating with a nearby bare mass m s to obtain the cor-responding η s mass, and then extracting the tuned mass using: m tuned0 s = m s m phys η s m η s ! . (A5)Our η s data are presented in Table V, which shows that thetuned mass is quite insensitive to small variations in m s . Wedo not have η s results for ensemble 7; there the tuned s massis based on an interpolation between results from ensemble 8and another ensemble that has similar parameters but with am ℓ = 0 . .Table II shows that m tuned0 c is more accurate than m tuned0 s .This is because the uncertainties in the value of the latticespacing have a smaller impact on the c mass because thecancellation described above only happens for heavy quarks(where m η h ≈ m h ).We set the u and d masses equal to their average, m ℓ ≡ m u + m d , (A6)and set m ℓ equal to the tuned s mass (above) divided by thephysical value of the quark mass ratio [33] m s m ℓ = 27 . . (A7) α MS ( µ, δm sea ) and a ( δm sea ) The beta function in the MS scheme is, by definition, inde-pendent of sea-quark masses. Thus the coupling’s evolution isunchanged by detuned sea-quark masses — dα MS ( µ, δm sea ) d log µ = β ( α MS ( µ, δm sea )) (A8) TABLE V. Simulation results for the η s mass am η s corresponding todifferent values of the bare s mass am s and different gluon ensem-bles. The ensembles are described in Table II, although we use manymore configurations for our η s analysis than are indicated there. Es-timates for the tuned bare s mass (Eq. (A5)) are also given.ensemble am s am η s am tuned0 s — but mass dependence enters through the low-energy start-ing point for that evolution implied by the scale-setting pro-cedure used in the lattice simulation. Such mass dependencecan enter only through an overall renormalization of the scaleparameter µ : α MS ( µ, δm sea ) = α MS ( ξ α µ ) (A9)where α MS ( µ ) ≡ α MS ( µ, δm sea = 0) (A10)is the MS coupling for physical sea-quark masses. The scalefactor, ξ α ≡ g α δm sea uds m s + g a ,α δm sea uds m s (cid:18) m c π/a (cid:19) + g c,α δm sea c m c + O ( δm ) , (A11)depends upon the differences between the masses m q used inthe simulation and the tuned values of those masses m tuned q (Table II and Sec. A 1): δm sea uds ≡ X q = u,d,s (cid:0) m q − m tuned q (cid:1) (A12) δm sea c ≡ m c − m tuned c . (A13)Function α MS ( ξ α µ ) satisfies the standard evolution equation(Eq. (A8)) because ξ α is independent of µ .We work to first order in δm sea because higher-order termsare negligible in our simulations. As suggested above, heleading-order dependence is particularly simple because weuse iso-singlet mesons ( η c and η s ) to set the c and s masses;in particular, there are no chiral logarithms of the u/d mass inleading order.4We expect coefficients g α and g a ,α in ξ α to be of or-der / since corrections linear in light-quark masses mustbe due to chiral symmetry breaking and so should be of or-der δm sea / Λ where Λ ≈ m s . As we discuss below, g c,α can be estimated from perturbation theory and is again of or-der / . We treat these coefficients as fit parameters in ouranalysis, with priors: g α = 0 ± . , g a ,α = 0 ± . , g c,α = 0 ± . . (A14)The rescaling factor ξ α is closely related to the dependenceof the lattice spacing on the sea-quark masses used in the sim-ulation. The lattice spacing is primarily a function of the barecoupling α lat used in the lattice action, but it also varies withthe sea-quark masses, in our scheme, when the bare couplingis held constant. As discussed above, this is because of sea-mass dependence in the quantity used to define the latticespacing, a/w in our case. The relationship with ξ α can beunderstood by examining the MS coupling at scale µ = π/a .There it is related to the bare coupling by a perturbative ex-pansion, α MS ( π/a, δm sea ) = α MS ( ξ α π/a )= α lat + ∞ X n =2 c MS n α n lat , (A15)that is mass-independent up to corrections of O (( am c ) α s ) ,which are negligible in our analysis. This formula implies that α MS ( ξ α π/a ) is constant if α lat is, and therefore that ξ α /a must be constant as well. Consequently the lattice spacingmust vary with δm sea like a ( δm sea ) ≈ ξ α a phys (A16)if the bare coupling is held constant, where a phys is the latticespacing when the sea-quark masses are tuned to their physicalvalues — that is, a phys ≡ a ( δm sea = 0) .We use this variation in the lattice spacing to read off theparameters in ξ α . Our simulation results fall into four groupsof gluon ensembles, with lattice spacings around 0.15 fm,0.12 fm, 0.09 fm and 0.06 fm. Each group corresponds to asingle value of the bare lattice coupling α lat , and several dif-ferent values of light sea-quark mass. Within a single group,then, the values we obtain for a/w from our simulationsshould vary as ( a/w ) sim = ξ α × ( a/w ) phys , (A17)where the parameters g α , g a ,α and g c,α in ξ α (Eq. (A11)) arethe same for all four groups of data.We fit our simulation results for a/w , simultaneously forall four groups, as functions of g α , g a ,α and g c,α . We alsotreat the value of ( a/w ) phys for each group as a fit parameter.The resulting fit is shown in Fig. 8 where we plot ( a/w ) sim ( a/w ) phys versus δm sea uds /m s . − . . . . . . δ m sea uds / m s . . . . . . ( a / w ) s i m / ( a / w ) phy s FIG. 8. The ratio of the simulation lattice spacing with detuned sea-quark masses to the lattice spacing with physical sea-quark massesas a function of the light-quark mass detuning (in units of the s quarkmass). Results are shown for four different sets of data, each corre-sponding to a different bare lattice coupling. The approximate lat-tice spacings for these sets are: 0.15 fm (red points), 0.12 fm (cyan),0.09 fm (green), and 0.06 fm (blue). The dashed line and gray bandshow the mean and standard deviation of our best fit to these data.The fit has a χ per degree of freedom of 0.23 for 9 degrees of free-dom ( p -value of 0.99). The fit is excellent, and shows that g α = 0 . . Our fitis not very sensitive to g a ,α and g c,α — their impact on ξ α istoo small — and gives results for these that are essentially thesame as the prior values. m h ( µ, δm sea ) and m c ( δm sea ) The evolution equations for the heavy quark’s MS mass areunchanged by sea-mass detunings: d log( m h ( µ, δm sea )) d log µ = γ m ( α MS ( µ, δm sea )) (A18)Consequently any sea-mass dependence must enter throughrescalings: m h ( µ, δm sea ) = ξ m m h ( ξ α µ ) (A19)where ξ α is defined above (Eq. (A11)), ξ m is independentof µ , and m h ( µ ) ≡ m h ( µ, δm sea = 0) (A20)is the MS mass for physical sea-quark masses. We parame-terize ξ m similarly to ξ α but allowing for the coefficients todepend upon the heavy-quark mass: ξ m = 1 + g m ( m η h /m η c ) ζ δm sea uds m s + g a ,m ( m η h /m η c ) ζ δm sea uds m s (cid:18) m c π/a (cid:19) + · · · (A21)5Again we expect g m and g a ,m to be of order / , and wetreat them as fit parameters with priors: g m = 0 ± . , g a ,m = 0 ± . . (A22)We parameterize the dependence on heavy-quark mass withthe factors ( m η h /m η c ) ζ where ζ is a fit parameter with prior: ζ = 0 ± . (A23)The sea-mass dependence in ξ m comes from the quantityused to tune the heavy-quark mass in simulations. We tunethese masses to give the correct physical mass for η h — thatis, the mass obtained when the sea-quark masses are tunedto their physical values and the lattice spacing is set to zero.This means that any sea-mass dependence in m η h is pushedinto the rescaling factor ξ m in Eq. (A19). The physical sizeof η h mesons decreases as m η h increases, and this decreasesthe coupling with light sea-quarks. Thus we expect ζ > in Eq. (A21); our fit finds ζ = 0 . .In principle, ξ m should depend upon δm sea c , as wellas δm sea uds . Perturbation theory, however, indicates that thisdependence is negligible in our simulations. Thus we haveomitted such terms from ξ m . We have verified that they arenegligible by comparing fits that include δm sea c terms with thefit without them.The rescaling factor ξ m is closely related to the sea-massdependence of the heavy quark’s bare mass, in much the sameway ξ α is related to the lattice spacing. The bare mass m h isproportional to the MS mass evaluated at µ = π/a : m h ∝ m h ( π/a, δm sea ) ∝ ξ m m h ( ξ α π/a ) . (A24)Since ξ α /a is sea-mass independent, we see that m h is pro-portional to ξ m , m h ( δm sea ) = ξ m m phys0 h , (A25)when the sea-quark masses are varied while holding the barecoupling fixed.This variation can be used to determine the parametersin ξ m , again in analogy to the previous section. As discussedin the previous section, our ensembles fall into four groupseach corresponding to a different value of the bare couplingconstant α lat . The masses am tuned0 c for each ensemble in Ta-ble II are tuned to give the physical η c mass for that ensemble.Therefore, within each group of ensembles, we expect am tuned0 c = ξ α ξ m × ( am c ) phys (A26)where ( am c ) phys is the value for properly tuned sea-quarkmasses.We fit our simulation results for am tuned0 c as functions of g m , g a ,m , g α , g a ,α , and g c,α . We use best-fit values from thefit in the previous section as priors for the last three of these fitparameters. The values of ( am c ) phys for the different groupsof ensembles are also fit parameters.The resulting fit is shown in Fig. 9, where we plot am tuned0 c / ( am c ) phys as a function of δm sea uds /m s . The fit isexcellent and shows that g m = 0 . , while g a ,m is es-sentially unchanged from its prior value (because our data arenot sufficiently accurate). − . . . . . . δ m sea uds / m s . . . . . . a m t un e d0 c / ( a m c ) phy s FIG. 9. The ratio of the bare c mass in lattice units used in the simu-lations to the bare mass with physical sea-quark masses as a functionof the light-quark mass detuning (in units of the s quark mass). Re-sults are shown for four different sets of data, each corresponding toa different bare lattice coupling. The approximate lattice spacings forthese sets are: 0.15 fm (red points), 0.12 fm (cyan), 0.09 fm (green),and 0.06 fm (blue). The dashed line and gray band show the meanand standard deviation of our best fit to these data. The fit has a χ per degree of freedom of 0.15 for 9 degrees of freedom ( p -value of1.0). c Quarks and Decoupling
Heavy quarks decouple from low-energy physics, andtherefore variations in δm sea c should have no impact onphysics (like w ) that probes momentum scales smallerthan m c . We can, however, introduce (apparent) violationsof the decoupling theorem through the scheme used to set thelattice spacing. In particular, decoupling is violated by anyscheme that holds the lattice spacing fixed (together with thebare coupling α lat ) as δm sea c is varied. On the contrary, decou-pling is preserved by schemes that hold a low-energy ( < m c )quantity like w fixed, instead of the lattice spacing [56].The difference between these schemes arises because therunning of the QCD coupling is modified in a detuned theoryfor scales between m sea c and m sea c + δm sea c , resulting in a mis-match between low and high energy values of the coupling.Physics below m c is determined by the n f = 3 coupling con-stant, which, by decoupling, should be independent of δm sea c .To see how this works, we examine lowest-order perturba-tion theory where α ( n f ) s ( µ ) = 2 πβ ( n f ) log( µ/ Λ ( n f ) ) (A27)with β ( n f ) ≡ − n f / , and α (3) s ( µ ) = α (4) s ( µ, δm sea c ) (A28)at µ = m c + δm sea c . Here Λ (3) must be independent of δm sea c ,by decoupling, while Λ (4) must vary with δm sea c to cancel theeffect of the shift in the match point µ = m c + δm sea c . It is6straightforward to show that Λ (4) ( δm sea c ) ≈ m c (cid:18) Λ (3) m c (cid:19) β (3) /β (4) (cid:18) − δm sea c m c (cid:19) ≈ Λ (4)phys × (cid:18) − δm sea c m c (cid:19) (A29)where Λ (4)phys is the value for physical sea-quark masses. Thusthe decoupling theorem requires that α (4) s ( µ, δm sea c ) = α (4) s (cid:18) µ × (cid:18) δm sea c m c (cid:19)(cid:19) . (A30)By comparing with Eqs. (A9) and (A11), we see that g c,α = 225 + O ( α s ) , (A31)and, therefore, that the lattice spacing varies with δm sea c (Eq. (A16)).There is an analogous effect in the heavy-quark mass, butthe mass dependence in ξ m is suppressed by α s and so is neg-ligible in our analysis.This analysis shows that a constant lattice spacing is in-compatible with the decoupling theorem. The scheme we useavoids this problem by allowing the lattice spacing to varywith δm sea c , while holding the value of w constant (as re-quired by the decoupling theorem applied to w itself). Theviolation of the decoupling theorem in the former case is onlyapparent; results from all schemes should agree when the sea-quark masses are tuned to their physical values. Appendix B: Previous Method
The analysis in our previous ( n f = 3 ) paper used a differentdefinition for the reduced moments with n ≥ : R n ≥ = m η h m h (cid:16) G n /G (0) n (cid:17) / ( n − (B1)instead of Eq. (3). As a result these moments equal z ( m η h , µ ) r n ( α MS , µ ) in perturbation theory where z ( m η c , µ ) ≡ m η h m h ( µ ) (B2) replaces z c ( µ ) , which is defined at the c mass instead of m h .Fits to these moments give both the coupling and the function z ( m η h , µ ) , from which the c and b masses can be extracted.We analyzed our data using the old definition, parameter-izing the m η h dependence of z ( m η c , µ ) with a cubic spline.The values for the R n moments used are given in Table VI.We obtained results that agree with the results obtained fromour new method to within a standard deviation, but are notquite as accurate: α MS (5 GeV , n f = 4) = 0 . (B3) m c (3 GeV , n f = 4) = 0 . . (B4) TABLE VI. Simulations results for η h masses and reduced moments R n (old definition) with various bare heavy-quark masses am h andgluon ensembles (first column, see Table II). Data from gluon en-sembles 1–3 are not listed because they were not used in the analysisin Appendix B. am h am η h R R R R The older method is more complicated because it attemptsto determine the coupling at the same time as it determinesthe functional dependence of z ( m η h , µ = 3 m h ) . In the newmethod, z ( m η h , µ = 3 m h ) is replaced by z c ( µ ) , whose de-pendence on µ is known a priori from perturbative QCD. [1] For a review see G. P. Lepage, P. B. Mackenzie, and M. E.Peskin, (2014), arXiv:1404.0319 [hep-ph].[2] C. McNeile, C. Davies, E. Follana, K. Hornbostel, and G. Lep-age (HPQCD Collaboration), Phys.Rev. D82 , 034512 (2010),arXiv:1004.4285 [hep-lat].[3] E. Follana et al. (HPQCD Collaboration, UKQCDCollaboration), Phys.Rev.
D75 , 054502 (2007), arXiv:hep-lat/0610092 [hep-lat].[4] A. Hart, G. von Hippel, and R. Horgan (HPQCDCollaboration), Phys.Rev.
D79 , 074008 (2009),arXiv:0812.0503 [hep-lat].[5] M. A. Shifman, A. Vainshtein, and V. I. Zakharov,Nucl.Phys.
B147 , 385 (1979).[6] K. Chetyrkin, J. H. Kuhn, and C. Sturm, Eur.Phys.J.
C48 , 107 (2006), arXiv:hep-ph/0604234 [hep-ph].[7] R. Boughezal, M. Czakon, andT. Schutzmeier, Phys.Rev.
D74 , 074006 (2006),arXiv:hep-ph/0605023 [hep-ph].[8] A. Maier, P. Maierhofer, and P. Marqaurd,Phys.Lett.
B669 , 88 (2008), arXiv:0806.3405 [hep-ph].[9] Y. Kiyo, A. Maier, P. Maierhofer, and P. Marquard,Nucl.Phys.
B823 , 269 (2009), arXiv:0907.2120 [hep-ph].[10] A. Maier, P. Maierhofer, P. Marquard, and A. Smirnov,Nucl.Phys.
B824 , 1 (2010), arXiv:0907.2117 [hep-ph].[11] D. J. Broadhurst, P. Baikov, V. Ilyin, J. Fleis-cher, O. Tarasov, et al. , Phys.Lett.
B329 , 103 (1994),arXiv:hep-ph/9403274 [hep-ph].[12] G. Martinelli, C. Pittori, C. T. Sachrajda, M. Testa,and A. Vladikas, Nucl.Phys.
B445 , 81 (1995),arXiv:hep-lat/9411010 [hep-lat].[13] A. Bazavov et al. (MILC collaboration),Phys.Rev.
D82 , 074501 (2010), arXiv:1004.0342 [hep-lat].[14] A. Bazavov et al. (MILC Collaboration),Phys.Rev.
D87 , 054505 (2013), arXiv:1212.4768 [hep-lat].[15] G. Lepage, B. Clark, C. Davies, K. Hornbostel,P. Mackenzie, et al. , Nucl.Phys.Proc.Suppl. , 12 (2002),arXiv:hep-lat/0110175 [hep-lat].[16] S. Borsanyi, S. Durr, Z. Fodor, C. Hoelbling, S. D. Katz, et al. ,JHEP , 010 (2012), arXiv:1203.4469 [hep-lat].[17] R. Dowdall, C. Davies, G. Lepage, and C. Mc-Neile (HPQCD Collaboration), Phys.Rev.
D88 , 074504 (2013),arXiv:1303.1670 [hep-lat].[18] T. van Ritbergen, J. Vermaseren, and S. Larin,Phys.Lett.
B400 , 379 (1997), arXiv:hep-ph/9701390 [hep-ph].[19] M. Czakon, Nucl.Phys.
B710 , 485 (2005),arXiv:hep-ph/0411261 [hep-ph].[20] J. Vermaseren, S. Larin, and T. van Ritbergen,Phys.Lett.
B405 , 327 (1997), arXiv:hep-ph/9703284 [hep-ph].[21] K. Chetyrkin, Phys.Lett.
B404 , 161 (1997),arXiv:hep-ph/9703278 [hep-ph].[22] V. Novikov, M. A. Shifman, A. Vainshtein, and V. I. Zakharov,Nucl.Phys.
B249 , 445 (1985).[23] M. A. Shifman, Prog.Theor.Phys.Suppl. , 1 (1998),arXiv:hep-ph/9802214 [hep-ph].[24] See, for example, M. Shifman, (2013),arXiv:1310.1966 [hep-th].[25] S. Hollands and C. Kopper,Commun.Math.Phys. , 257 (2012),arXiv:1105.3375 [hep-th].[26] D. Pappadopulo, S. Rychkov, J. Espin, and R. Rattazzi,Phys.Rev.
D86 , 105043 (2012), arXiv:1208.6449 [hep-th].[27] K. Chetyrkin, B. A. Kniehl, and M. Steinhauser,Nucl.Phys.
B510 , 61 (1998), arXiv:hep-ph/9708255 [hep-ph].[28] K. Olive et al. (Particle Data Group),Chin.Phys.
C38 , 090001 (2014) and 2013 partial updatefor the 2014 edition.[29] C. McNeile, A. Bazavov, C. Davies, R. Dowdall,K. Hornbostel, et al. , Phys.Rev.
D87 , 034503 (2013),arXiv:1211.6577 [hep-lat].[30] I. Allison et al. , Phys.Rev.
D78 , 054513 (2008),arXiv:0805.2999 [hep-lat].[31] The precise definition of our error budgets is described in Ap-pendix A of C. Bouchard, G. P. Lepage, C. Monahan, H. Na,and J. Shigemitsu, (2014), arXiv:1406.2279 [hep-lat].[32] C. Davies, C. McNeile, K. Wong, E. Follana, R. Horgan, et al. (HPQCD Collaboration), Phys.Rev.Lett. , 132003 (2010),arXiv:0910.3102 [hep-ph]. [33] A. Bazavov et al. (Fermilab Lattice and MILC Collaborations),(2014), arXiv:1407.3772 [hep-lat].[34] K. Chetyrkin, J. Kuhn, A. Maier, P. Maierhofer,P. Marquard, et al. , Phys.Rev.
D80 , 074010 (2009),arXiv:0907.2110 [hep-ph].[35] K. Chetyrkin, J. Kuhn, A. Maier, P. Maierhofer,P. Marquard, et al. , Theor.Math.Phys. , 217 (2012),arXiv:1010.6157 [hep-ph].[36] N. Carrasco, A. Deuzeman, P. Dimopoulos, R. Frezzotti,V. Gimenez, et al. , (2014), arXiv:1403.4504 [hep-lat].[37] S. Durr and G. Koutsou, Phys.Rev.Lett. , 122003 (2012),arXiv:1108.1650 [hep-lat].[38] B. Blossier et al. (ETM Collaboration),Phys.Rev.
D82 , 114513 (2010), arXiv:1010.3659 [hep-lat].[39] S. Durr, Z. Fodor, C. Hoelbling, S. Katz, S. Krieg, et al. ,Phys.Lett.
B701 , 265 (2011), arXiv:1011.2403 [hep-lat].[40] R. Arthur et al. (RBC Collaboration, UKQCD Collaboration),Phys.Rev.
D87 , 094514 (2013), arXiv:1208.4412 [hep-lat].[41] H. Georgi and C. Jarlskog, Phys.Lett.
B86 , 297 (1979).[42] C. Davies et al. (HPQCD Collaboration),Phys.Rev.
D78 , 114507 (2008), arXiv:0807.1687 [hep-lat].[43] E. Shintani, S. Aoki, H. Fukaya, S. Hashimoto,T. Kaneko, et al. , Phys.Rev.
D82 , 074505 (2010),arXiv:1002.0371 [hep-lat].[44] A. Bazavov, N. Brambilla, X. Garcia i Tormo, P. Pe-treczky, J. Soto, et al. , Phys.Rev.
D86 , 114031 (2012),arXiv:1205.6155 [hep-ph].[45] K. Jansen, F. Karbstein, A. Nagy, and M. Wag-ner (ETM Collaboration), JHEP , 025 (2012),arXiv:1110.6859 [hep-ph].[46] P. Fritzsch, F. Knechtli, B. Leder, M. Marinkovic,S. Schaefer, et al. , Nucl.Phys.
B865 , 397 (2012),arXiv:1205.5380 [hep-lat].[47] B. Blossier et al. (ETM Collaboration),Phys.Rev.
D89 , 014507 (2014), arXiv:1310.3763 [hep-ph].[48] A. Lee et al. (HPQCD Collaboration),Phys.Rev.
D87 , 074018 (2013), arXiv:1302.3739 [hep-lat].[49] A. X. El-Khadra, A. S. Kronfeld, and P. B. Mackenzie,Phys.Rev.
D55 , 3933 (1997), arXiv:hep-lat/9604004 [hep-lat].[50] B. Colquhoun, R. Dowdall, C. Davies, K. Horn-bostel, and G. Lepage (HPQCD Collaboration), (2014),arXiv:1408.5768 [hep-lat].[51] Perturbation theory suggests that the important scales in w are of order Q w = 1 / p w ≈ MeV. See M. Luscher,JHEP , 071 (2010), arXiv:1006.4518 [hep-lat].[52] O. B¨ar and M. Golterman, Phys.Rev.
D89 , 034505 (2014),arXiv:1312.4999 [hep-lat].[53] C. Davies, C. McNeile, E. Follana, G. Lepage, H. Na, et al. (HPQCD Collaboration), Phys.Rev.
D82 , 114504 (2010),arXiv:1008.4018 [hep-lat].[54] E. B. Gregory, C. T. Davies, I. D. Kendall, J. Ko-ponen, K. Wong, et al. (HPQCD Collaboration),Phys.Rev.
D83 , 014506 (2011), arXiv:1010.3848 [hep-lat].[55] C. Davies, K. Hornbostel, A. Langnau, G. Lep-age, A. Lidsey, et al. , Phys.Rev.Lett. , 2654 (1994),arXiv:hep-lat/9404012 [hep-lat].[56] Dimensionless ratios of low-energy quantities are indepen-dent of the lattice-spacing scheme, and must be independentof m c by the decoupling theorem. This means that a schemethat makes any one low-energy quantity — for example, w —independent of m c makes all other low-energy quantities inde-pendent of m cc