Determination of the charm quark mass in lattice QCD with 2+1 flavours on fine lattices
MMS-TP-21-01
Determination of the charm quark mass in lattice QCDwith flavours on fine lattices
LPHA A Collaboration
Jochen Heitger , Fabian Joswig and Simon Kuberski † Westfälische Wilhelms-Universität Münster, Institut für Theoretische Physik,Wilhelm-Klemm-Straße 9, 48149 Münster, Germany
Abstract
We present a determination of the charm quark mass in lattice QCD with three activequark flavours. The calculation is based on PCAC masses extracted from N f = 2 + 1flavour gauge field ensembles at five different lattice spacings in a range from 0.087 fm downto 0.039 fm. The lattice action consists of the O( a ) improved Wilson-clover action and atree-level improved Symanzik gauge action. Quark masses are non-perturbatively O( a )improved employing the Symanzik-counterterms available for this discretisation of QCD.To relate the bare mass at a specified low-energy scale with the renormalisation groupinvariant mass in the continuum limit, we use the non-pertubatively known factors thataccount for the running of the quark masses as well as for their renormalisation at hadronicscales. We obtain the renormalisation group invariant charm quark mass at the physicalpoint of the three-flavour theory to be M c = 1486(21) MeV. Combining this result withfive-loop perturbation theory and the corresponding decoupling relations in the MS scheme,one arrives at a result for the renormalisation group invariant charm quark mass in thefour-flavour theory of M c ( N f = 4) = 1548(23) MeV. In the MS scheme, and at finite energyscales conventional in phenomenology, we quote m MSc ( m MSc ; N f = 4) = 1296(19) MeV and m MSc (3 GeV; N f = 4) = 1007(16) MeV for the renormalised charm quark mass. Keywords:
Lattice QCD, Non-perturbative effects, Quark masses, Charm † Present address:
Helmholtz-Institut Mainz, Johannes Gutenberg-Universität Mainz, 55099 Mainz,Germany a r X i v : . [ h e p - l a t ] J a n Introduction
Apart from the strong coupling constant, quark masses are the only fundamental parametersof quantum chromodynamics (QCD), the theory of strong interaction. As such, their preciseknowledge is not only of general interest, but also essential in the search for new physicsbeyond the Standard Model of particle physics. In particular, heavy quark masses serve askey parametric inputs for tests of the Standard Model, via hadronic contributions by QCDmatrix elements to weak decays and ensuing constraints on CKM matrix elements [1], aswell as in Higgs branching ratios to charm and bottom quarks; see for instance Refs. [2, 3].The extraction of these and other Standard Model parameters that have their origin in thestrong interaction requires to evaluate the QCD path integral. Since quarks are confinedin bound states and cannot propagate freely, their properties have to be investigatedin the hadronic regime at low energies where perturbative methods fail. We base ourcalculations on lattice QCD, where the theory is formulated on a discretised space-timegrid to enable, in conjunction with Monte Carlo simulation techniques, an ab initio andthereby non-perturbative solution of the underlying highly non-linear equations with, inprinciple, full control of statistical uncertainties and systematic effects.Our determination of the charm quark mass is performed on a subset of the gaugeconfiguration ensembles, generated by the
Coordinated Lattice Simulations (CLS) coopera-tive effort of European lattice QCD teams, with N f = 2 + 1 flavours of non-perturbativelyO( a ) improved Wilson quarks and the tree-level Symanzik-improved gauge action [4, 5].The sea quark content in these discretised QCD configurations comprises a doublet of lightdegenerate quarks (with masses to be identified with the mean up and down quark mass)and a third, single flavour representing the strange quark. The masses employed in ournumerical calculations differ from their physical counterparts in nature, which requiresapproaching the physical point through a chiral extrapolation to the physical quark massvalues corresponding to QCD with 2 + 1 flavours.In practice, the tuning of the parameters of the theory is guided by the computationof meson masses, where the physical point is defined by the pion and kaon (as well as somecharmed meson) masses taking their physical values. As the sum of sea quark masses isheld constant in our setup, the strange quark mass reaches its physical value from below.However, we would like to anticipate already now that no significant dependence on thelight quark masses is observed in our analysis of observables in the charm sector understudy.The charm quark is included as partially quenched valence quark in our calculation,which means that it is not incorporated as a dynamical degree of freedom in the latticeaction governing the path integral of the theory. The size of the systematic uncertaintyinduced by this approximation has been examined in Ref. [6] and was found to be, owingto the decoupling of heavier quarks, of the order of a few per cent at most. On the otherhand, the comparison of lattice QCD charm quark determinations relying on 2 + 1 [7–9]and 2 + 1 + 1 dynamical flavours [10–14], collected and reviewed in Ref. [15], demonstratesan excellent agreement. As for the effects of QED and isospin breaking, the statisticalerror on our final result will turn out to be substantially larger than the level, at whichone expects their influence to affect current lattice QCD computations; this justifies toneglect them here.Lattice QCD offers a sound, first principles definition of quark masses that derives2rom a lattice version of the low-energy theorem of QCD known as the PCAC, viz., partiallyconserved axial current relation, holding in the chiral symmetry realm of the theory. Oncebare quark masses based on the PCAC relation on the lattice are extracted (as we do it herein case of the charm quark mass), they need to be properly renormalised for their continuumlimit to exist. A convenient object to consider is the renormalisation group invariant (RGI)quark mass, as it is a renormalisation scheme- and scale-independent quantity (if therenormalisation conditions are imposed at zero quark mass [16]), and the determination ofthe RGI charm quark mass is thus the main goal of this paper. For connecting the bare(PCAC) mass at a specified low-energy scale with the RGI mass in the continuum limit, wecan draw on the study exposed in Ref. [17], where the renormalisation group running of thequark mass over a wide range of scales for the present three-flavour discretisation of QCDwas performed with full non–perturbative precision in the Schrödinger functional scheme,following the strategy of Ref. [18]. The results of Ref. [17] split into a factor accounting forthe universal (i.e., regularisation and flavour independent) running of the quark massesand an associated renormalisation factor at a low-energy scale in an appropriate hadronicscheme, defined via the Schrödinger functional. Moreover, our computations involvenon-perturbatively determined improvement coefficients and renormalisation constantsfrom Refs. [17, 19–22], implementing Symanzik improvement for removal of discretisationeffects from correlation functions and chiral Ward identities for finite normalisations ofquark bilinears (such as the axial vector current) in the O( a ) theory. Therefore, leadingdiscretisation errors of O( a ) in the bulk (except for O( g a ) uncertainties at the timeboundaries) are encountered when taking the continuum limit of correlation functions andrenormalised quark masses deduced from them.The combination of non-perturbative O( a ) improvement and the range of latticespacings from ≈ .
087 fm down to very fine resolutions of ≈ .
039 fm covered in our work(where the calibration of the lattice spacing to its physical value was done through theintermediate hadronic gradient flow scale t in Ref. [23]) allow us to carry out a controlledcontinuum limit extrapolation, together with the chiral extrapolation in a joint global fitalong a trajectory of constant physics towards the physical point of the theory. In thiscontext, the use of different lattice definitions of the charm quark mass, which distinguishthemselves considerably in their cutoff effects, proves to be very valuable, because we areable to confirm that no significant systematic uncertainties are introduced by the globalextrapolation to the joint chiral and continuum limit. We find that linear cutoff effects areabsent (as expected), but those of O( a ) and higher give relevant contributions already atlattice spacings around 0 .
06 fm.The structure of this paper is as follows. In Section 2, we introduce the three differentdefinitions of the bare quark mass on the lattice employed and explain their renormalisationpattern in the O( a ) improved theory. Section 3 specifies the set of gauge field configurationensembles and how they realise, in lattice parameter space, the approach to the point ofphysical quark masses (fixed, in practice, by physical values of pion, kaon and suitablecharmed meson masses) along a chiral trajectory and in the continuum limit. There, wealso detail the computation of two-point correlation functions, as well as the extractionof meson and bare quark masses from them, and provide the expressions from which therenormalised charm quark mass is obtained and eventually translated into the associatedRGI mass. The joint global chiral-continuum extrapolation of the resulting estimates atstill unphysical meson masses and finite lattice spacings are then thoroughly described3n Section 4, where we also elucidate the model averaging procedure that is adopted toarrive at our final result and to quantify its systematic uncertainties. Our final result onthe RGI charm quark mass is presented in Section 5, together with a careful discussion ofthe various contributions from the single steps of our calculation to its full error budget.In Section 6, we conclude by summarising our findings and compare our final value forthe charm quark mass with existing ones from other lattice QCD determinations, seealso Appendix A, where the conversion through perturbative five-loop running to otherconventional energy scales in the MS renormalisation scheme is outlined, which representsthe customary reference scheme in phenomenological applications. Lastly, results fromseveral intermediate steps of the analysis, such as meson and (bare) quark masses in latticeunits, are tabulated in Appendix B.A preliminary account of our work has been given in Ref. [24]. The determinationof the mass-degenerate light and strange quark masses on (2 + 1)-flavour CLS ensembles,along the same renormalisation strategy and partly similar analysis methods, was presentedin Ref. [25]. Our determination of the physical charm quark mass is based on various different definitionsof the quark mass in the formulation of lattice QCD with Wilson quarks. We begin ourdiscussion of these definitions with the bare subtracted quark mass of flavour i , which isdefined as m q ,i = m ,i − m cr ≡ a (cid:18) κ i − κ cr (cid:19) , (2.1)where κ i is the corresponding hopping parameter and κ cr is the critical hopping parameterdefined from vanishing current quark masses in the SU ( N f ) symmetric limit. The additivequark mass renormalisation via m cr is required due to the loss of chiral symmetry in theWilson fermion action. The average of two subtracted quark masses of flavours i and j is m q ,ij ≡
12 ( m q ,i + m q ,j ) . (2.2)Since composite operators built from Wilson fermion fields show considerable cutoffeffects that are enhanced at the scale of the charm quark, O( a ) improvement à la Symanzikdrastically improves their approach towards the continuum limit. A detailed discussion ofthe underlying improvement pattern in the case of three flavours of quarks can be found inRef. [26]. We repeat the expression for the renormalised and improved subtracted quarkmass, m i, R ≡ Z m (cid:26)(cid:20) m q ,i + ( r m −
1) Tr [ M q ] N f (cid:21) + aB i (cid:27) , (2.3) B i = b m m ,i + ¯ b m m q ,i Tr [ M q ]+ ( r m d m − b m ) Tr (cid:2) M (cid:3) N f (2.4)+ ( r m ¯ d m − ¯ b m ) (cid:0) Tr [ M q ] (cid:1) N f , M q = diag( m q , , m q , , . . . , m q ,N f ) , (2.5)and the mass renormalisation constant Z m that not only depends on the bare coupling, g , but also explicitly on the renormalisation scheme and associated renormalisation scale aµ . As apparent from eq. (2.3), the subtracted quark mass receives a shift from non-zerosea quark masses which scales with ( r m − r m is the ratio of flavour singlet andnon-singlet scalar density renormalisation constants. The cutoff effects at order O( a ) aremass dependent and get fully cleared if the coefficients of the corresponding terms are fixednon-perturbatively.As an alternative definition of the quark mass on the lattice we employ the barecurrent quark mass entering the PCAC (partially conserved axial current) relation. Fortwo distinct flavours i = j it reads m ij = h ( ˜ ∂ A ij ( x ) + ac A ∂ ∗ ∂ P ij ( x )) P ji ( y ) i h P ij ( x ) P ji ( y ) i , (2.6)in terms of the improved axial vector current( A I ) ijµ = A ijµ + ac A ∂ µ P ij , (2.7)and the pseudoscalar density P ij , which in the presence of non-degenerate masses areoff-diagonal bilinear fields. The standard choice for the discretised derivatives is given bythe forward and backward differences a∂ µ f ( x ) ≡ f ( x + a ˆ µ ) − f ( x ) , a∂ ∗ µ f ( x ) ≡ f ( x ) − f ( x − a ˆ µ ) , (2.8)while ˜ ∂ denotes the average of these two. Following Refs. [21, 27, 28], we also adopt anotherrealisation of discretised derivatives based on the substitutions˜ ∂ µ → ˜ ∂ µ (cid:16) − a ∂ ∗ µ ∂ µ (cid:17) , ∂ ∗ µ ∂ µ → ∂ ∗ µ ∂ µ (cid:16) − a ∂ ∗ µ ∂ µ (cid:17) . (2.9)Since the discretisation effects introduced by these derivatives acting on smooth functionsare of O( a ), we will refer to this second choice as improved derivatives.Upon renormalisation, eq. (2.6) leads to the renormalised O( a ) improved current quarkmass m R ,ij , m R ,ij = Z A Z P m ij h b A − b P ) am q ,ij + (¯ b A − ¯ b P ) a Tr [ M q ] i , (2.10)which in the theory with mass non-degenerate quarks is equivalent to the average of tworenormalised current quark masses of flavour i and j . Here, Z A ( g ) and Z P ( g , aµ ) labelthe renormalisation constants of the axial current and the pseudoscalar density, respectively,the latter carrying the renormalisation scheme and scale dependences of the quark mass asdoes Z m in eq. (2.3) above.We can now equate both expressions for the renormalised quark mass, eqs. (2.3) and(2.10), to derive m q ,ij = m ij Z − ( r m −
1) Tr [ M q ] N f + O( am ij , a Tr [ M q ]) , (2.11) For ease of notation, we refrain from explicitly writing out corrections of O( a ) or higher. Z is the scale independent combination of renormalisation constants Z ( g ) ≡ Z m ( g , aµ ) Z P ( g , aµ ) Z A ( g ) , (2.12)and employ this relation to eliminate the bare subtracted quark mass in eq. (2.10) in favourof PCAC masses to obtain m R ,ij = Z A Z P m ij " b A − ˜ b P ) am ij (2.13)+ (¯ b A − ¯ b P ) Zr m − (˜ b A − ˜ b P ) ( r m − r m N f ! aM sum , with ˜ b A − ˜ b P ≡ ( b A − b P ) Z , (2.14) M sum ≡ m + m + · · · + m ( N f − N f + m N f = Zr m Tr [ M q ] + O( a Tr [ M q ]) . (2.15)Note that in eq. (2.13) all dependence on m cr has been removed.Yet another definition of the renormalised quark mass with different cutoff effects maybe constructed along the so-called ratio-difference method [29], where we again combinethe two definitions for renormalised quark masses, eqs. (2.3) and (2.10), but in a differentway than before. Namely, we introduce a ratio of current quark masses, r , and a differenceof subtracted quark masses, d , through r ij ≡ m ii m jj , d ij ≡ am q ,i − am q ,j , (2.16)where i and i label two distinct but mass degenerate quark flavours. This setup is chosensuch that the multiplicative renormalisation of the current quark masses cancels in theratio, while the additive renormalisation of the bare subtracted quark masses (via m cr )and the leading-order dependence on Tr [ M q ] drop out in the difference. As explained indetail in Ref. [30], the O( a ) improved versions of r ij and d ij can be written as r I ,ij = r ij [1 + ( b A − b P ) d ij ] , (2.17) d I ,ij = d ij " b m d ij r ij + 1 r ij − b m d ij − r m ) M sum r m m jj ( r ij − N f + a ¯ b m M sum Zr m . (2.18)From the latter, an estimator for the renormalised quark mass is formed by m (rd) i, R = Z m r I ,ij d I ,ij r I ,ij − , (2.19)with the mass renormalisation factor Z m already introduced in eq. (2.3).For the improvement and renormalisation of the axial current we use the improvementcoefficient c A calculated in Ref. [19] and the values of Z A ( g ) determined within the chirallyrotated Schrödinger functional [20]. The scale dependent renormalisation constant Z P wascomputed in Ref. [17], together with the universal factor that relates quark masses at finite6enormalisation scale to their RGI counterparts. The renormalisation constant Z definedin eq. (2.12) was calculated for our setup in Ref. [21]. It can also be extracted exploitingthe numerical results on Z S /Z P presented in Ref. [30].The non-perturbative determination of the improvement coefficients ( b A − b P ) and b m has been performed in Ref. [21], along with Z , for two different lines of constant physics(LCP) realised in the Schrödinger functional framework for the range of lattice spacingsrelevant here. More precisely, there is a set of coefficients computed on an LCP adapted foruse in a fully massless setup (called ‘LCP-0’ in [21]), whereas, by contrast, a further set ofestimates (labelled ‘LCP-1’) for Z , ( b A − b P ) and b m refers to an LCP at fixed heavy quarkmass in the valence sector. Therefore, particularly the latter appears to be predestinedfor applications with heavy (valence) quarks, because mass dependent cutoff effects of allorders are expected to be accounted for in the improvement coefficients and renormalisationconstant from this LCP. In fact, the investigations in Refs. [31, 32] have shown that thescaling of meson observables involving heavy quarks towards the continuum limit benefitsfrom including these mass dependent cutoff in the LCP definition of the improvementcoefficients and renormalisation factor in question. This conjecture will be tested later,when we discuss the chiral-continuum extrapolation of our results.For the improvement coefficient (¯ b A − ¯ b P ) and ¯ b m multiplying cutoff effects dependingon the sea quark masses, non-perturbative results with errors of O(100%) were obtainedin Refs. [33, 34] using coordinate-space methods. Given the large uncertainties in thesedeterminations and the fact that the coefficients are of O( g ) in perturbation theory (asopposed to O( g ) for the coefficients b X mentioned before), we ignore the correspondingterms of O( a Tr [ M q ]) in our improvement procedure.Finally, non-perturbative results for r m are available only for a subset of the CLSensembles used in our study [35]. However, a non-perturbative determination over the fullrange of couplings is almost finished and will be published soon [22]. While we thus employthe results from this work here, we also note that the associated contributions from thesea quark masses are highly suppressed compared to the heavy valence quark mass effectsthat are removed non-perturbatively. The numerical study described in the following has been carried out on the (2 + 1)-flavourCLS gauge field configuration ensembles described in Refs. [4, 5]. Gluons are implementedvia the tree-level improved Lüscher-Weisz gauge action [36] and quarks are treated as cloverimproved Wilson fermions [37]. All ensembles employed in our computations feature openboundary conditions in time direction to prevent a freezing of the topological charge [38,39].Twisted mass reweighting [40] has been applied in the simulations for the doublet of lightsea quarks and the RHMC algorithm [41, 42] has been used for the single strange quark inthe sea.Charm quarks are introduced in the valence sector only, rendering our theory apartially quenched approximation of QCD. Whereas the effect of a dynamical charm quarkin low-energy purely gluonic (and light flavour) observables is expected to be suppressedbelow the available accuracy [43], this is different for observables at the scale of the charmquark, and therefore a systematic uncertainty is present owing to neglecting a dynamicalcharm in the sea in our final result for the physical charm quark mass. A conservative7pper bound of 5% for the relative size of this uncertainty was estimated in Refs. [6, 44].Ensembles with five different lattice spacings in a range from 0 .
087 fm down to 0 .
039 fm,together with non-perturbative O( a ) improvement, allow for a controlled continuum limit ex-trapolation. Moreover, these ensembles cover pion masses in the range [200 MeV ,
420 MeV]such that the dependence of our results for the charm quark mass on the light quark massesmay be carefully investigated. In general, however, the effect of unphysically large seaquark masses on observables at the scale of the charm quark is expected to be small. Weperform our calculations on a set of ensembles, for which the sum of the bare sea quarkmasses Tr [ M q ] is held fixed. In such as setup, the strange quark mass approaches itsphysical value from below when the light quark masses are lowered towards their physicalvalues.The strategy to fix an LCP trajectory in bare lattice parameter space within theTr [ M q ] = const . CLS ensembles and its exact specification are explained in great detailin Refs. [23, 45]. Since a fixed value of the bare trace of the quark mass matrix does notimply a constant renormalized trace of the quark mass matrix, the generated gauge fieldconfiguration ensembles do not lie strictly on an LCP. Therefore, and since the deviationsof Tr [ M R , q ] from a constant value have been observed to be larger than what would beexpected from O( a ) cutoff effects [23], the chiral trajectory has been redefined in terms ofa constant value of φ ≡ t (cid:18) m + 12 m π (cid:19) , (3.1)where the gluonic quantity t is defined from the Wilson flow [46] and m π and m K denotethe masses of the lightest mesons composed of light and light-strange quarks, respectively.To approach the physical point in a chiral extrapolation, the value of φ on each ensemblehas to be fixed to the physical value φ phys4 = 8 t phys0 (cid:20)(cid:16) m physK (cid:17) + 12 (cid:16) m phys π (cid:17) (cid:21) = 1 . , (3.2)based on the physical values of pion and kaon masses in isospin symmetric QCD as detailedin Ref. [47] and the physical value of t determined in Ref. [23]. In this appropriatelyadapted setup, the dependence of physical observables on the sea quark masses can nowbe parametrised via φ ≡ t m π , (3.3)while the physical point is recovered by an extrapolation to φ phys2 = 8 t phys0 (cid:16) m phys π (cid:17) = 0 . . (3.4)To ensure φ = φ phys4 on each ensemble, we resort to small quark mass shifts onto thechiral trajectory by means of a Taylor expansion of the entering correlation functions asdescribed in Ref. [23]. Following this reference, the derivative of an observable O withrespect to a change in the quark mass m q ,i of flavour i is given byd hOi d m q ,i = * ∂ O ∂m q ,i + − * O ∂S∂m q ,i + + hOi * ∂S∂m q ,i + . (3.5)8n order to carry out this shift at the level of the correlation functions employed in our work,their partial derivatives were numerically evaluated (simultaneously with the calculationof the correlators themselves) utilising the mesons program package [48], augmented bysets of measurements for the derivative terms of the action with respect to the quark massin (3.5) that we had to extend compared to what was already available as a result of [23].The shift of observables depending on N f = 2 + 1 quark masses is then accomplished bythe replacement hOi → hOi + N f X i =1 ∆ m q ,i d hOi d m q ,i , (3.6)and we choose the same value of ∆ m q ,i for all three quark flavours. In practice, the valueof ∆ m q ,i on each ensemble was obtained in an iterative procedure that stops when theshifted φ matches its physical value well below the statistical precision. As already indicated above, we have performed correlator measurements on 15 ensembleson the Tr [ M q ] = const . trajectory from the (2 + 1)-flavour CLS gauge configuration database, using the mesons code [48]. The specifications of these ensembles and their associatedstatistics in the number of molecular dynamic units (MDU) are summarised in Table 1.To extract ground state meson masses as well as quark masses through the PCAC relation,we consider zero-momentum two-point correlation functions defined by f rsOO ( x , y ) = − a L X x , y (cid:10) O rs ( x , x ) O ,sr ( y , y ) (cid:11) , (3.7)where r and s are flavour indices. y denotes the time coordinate of the source, i.e., thetime slice where the source has a non-vanishing norm, and x is the time coordinate of thesink. The operators O rs are quark bilinear covariant fields composed as O rs ( x ) = ¯ ψ r ( x ) Γ ψ s ( x ) , (3.8)in terms of Γ, representing the combination of Euclidean gamma matrices that read γ γ for the time component of the axial current, A , and γ for the pseudoscalar density, P . Since translational invariance in time direction is broken by the open boundary conditions,we fix the temporal source positions y to the two time slices at a and T − a as proposedin Ref. [49]. To achieve good statistical accuracy, we average over 16 U(1) noise sourcesper source position and exploit time reversal symmetry.The correlation functions f rsOO were evaluated for all combinations of light and strangeas well as for two choices of heavy valence quark flavours. For the latter, the hoppingparameters have been chosen such that they encompass the supposed hopping parameterof the physical charm quark κ c sufficiently close. An overview of these hopping parametersis given in Table 2. Distance preconditioning [50], in its variant as as sketched in Ref. [51](and implemented in the mesons measurement code [48]), has been employed for thecalculation of the heavy quark propagators, in order to ensure the stability and accuracyof the associated numerical inversions of the Dirac operator across the complete temporalextent of the lattice. While also γ k for the spatial components of the vector current V k was included in the set of evaluatedcorrelators, the vector meson mass does not enter our final analysis, see below. able 1: Overview of the gauge field configuration ensembles (labelled by ‘id’) used in this study.We list the geometry ( T × L ) and the bare parameters of the action as well as the resultingapproximate lattice spacings and meson masses. The statistics included is quantified by the numberof molecular dynamics units N MDU . We performed measurements every 4 MDU except for J303and J500, where the spacing between two measurements is 8 MDU. The exponential autocorrelationtime is defined from t as given by eq. (3.10). β a fm id La Ta κ l κ s m π MeV m K MeV N MDU τ exp MDU .
40 0 .
087 H101 32 96 0 . κ l
416 416 8064 19 . . . . . . . . . . .
46 0 .
076 H400 32 96 0 . κ l
421 421 4180 25 . .
55 0 .
064 H200 32 96 0 . κ l
418 418 8000 36 . . κ l
411 411 3596 36 . . . . . . . . . . .
70 0 .
050 N300 48 128 0 . κ l
418 418 8188 59 . . . . . . . .
85 0 .
039 J500 64 192 0 . κ l
409 409 6312 98 . . . . We follow two different strategies to tune to the physical charm quark mass. Firstly,since Tr [ M q ] is held constant in our framework, the flavour averaged meson mass m ¯D ≡ m D + 13 m D s , (3.9)built from D and D s mesons, is approximately constant as well. At the same time, thestatistical precision of the corresponding effective masses is convincingly good. Secondly,we consider the effective mass of the connected part of the η c meson, called m η c from nowon. The signal of this observable is clean, and the effect of neglecting the disconnectedpart is expected to be subleading [52], particularly so for singlet mesons composed ofheavy quarks. Hence we have performed all measurements including heavy quarks for bothcharmed hopping parameters and, to eventually obtain the physical charm quark mass,interpolate the meson mass results to their respective physical, i.e., experimental valuesjointly with a chiral-continuum extrapolation procedure, as will be explained in detail inthe next section. Let us mention that we also had tested the spin-flavour averaged mesonmass constructed from heavy-light and heavy-strange pseudoscalar and vector mesons [24]to fix the charm quark mass. However, even though one is able to extract the ground Since the spin-dependent first-order correction of HQET to the meson mass is cancelled in thisaverage [53, 54], one might a priori expect reduced cutoff effects from the use of this average. able 2: Hopping parameters of the partially quenched heavy valence quarks in the charm regionserving to implicitly fix κ c , provided by [55] in the context of determining leptonic decay constantsof D and D s mesons [51, 56]. β κ h κ h .
40 0 . . .
46 0 . . .
55 0 . . .
70 0 . . .
85 0 . . states of the vector mesons (despite the correlators being considerably noisier than in thepseudoscalar channel) on all ensembles, we did not see any improvement in the approachof our results to the continuum limit compared to the other two methods. Moreover, dueto the still relatively large statistical uncertainties of the vector meson masses, there isalso no gain in overall statistical precision such that we did not include the spin-flavouraveraged meson mass in our final analysis.We carry out our statistical error analysis using the Γ-method [57] deriving fromautocorrelation functions, supplemented by automatic differentiation for the error propaga-tion [58], and take into account the effect of the exponential tails in the autocorrelationfunctions owing to slow modes [59]. As estimate for the required exponential autocorrelationtime on the CLS ensembles we recourse to the formula τ exp = 14(3) t a (3.10)quoted in Ref. [4] for our set of ensembles. Ground state masses of pseudoscalar mesons are extracted from the effective mass am eff ( x ) = ln (cid:18) f PP ( x ) f PP ( x + a ) (cid:19) , (3.11)which forms a plateau for sufficiently large source-sink separations. The computation of pionand kaon masses on lattices with open boundaries has been discussed in Refs. [4,23,49], andwe employ the same methods to determine the plateau region, i.e., fits to the exponentialcorrections originating from boundary effects and short source-sink separations. Includingthe dominant contributions amounts to model the effective masses with an ansatz am eff ( x ) = am PS (cid:16) c e − E x + c e − E ( T − x ) (cid:17) , (3.12)where m PS gives the mass of the pseudoscalar meson in question, E = m − m PS thedifference between ground and first excited state energies, and E (cid:117) m PS . The energiesand the parameters c i are chosen to be free fit parameters. The same procedure is appliedto the effective mass of the η c meson. For heavy-light and heavy-strange mesons we restrictourselves to a region where the statistical fluctuations are small such that effects of the11oundary at maximal source-sink separation do not affect the quality of the plateau. Weshow representative effective masses for these mesons, entering into the interpolations tofix the physical charm quark mass, in the left part of Figure 1; there the exponentialcorrections for small source-sink separations and close to the boundaries are clearly visible.To specify the optimal fit window, we vary the minimal source-sink distance includedin the fit and monitor its quality through the so-called ‘corrected’ χ , cf. eq. (4.8), of anuncorrelated fit. χ is based on the correlations present in the data, as motivated andelaborated in Ref. [60]. The beginning of the plateau range is set to the time slice wherethe contribution of the excited state corrections is four times smaller than the statisticalerror on the effective mass. The upper end of the fit interval is taken from a fit to theboundary corrections in case of mesons with constant signal to noise ratio, while it is abound of 3% on the relative error of the correlation function in case of heavy-light mesons.The meson mass can then be defined either as the weighted plateau average within theplateau range thus determined or directly from the fit to the plateau and its exponentialcorrections according to eq. (3.12). Because of their compatibility with each other, inthe final analysis we just opted for the latter. Note that thanks to the use of distancepreconditioning the data on heavy quark propagators over the full temporal range can beincorporated in the plateau fits. This effect is actually enhanced on the J ensembles closeto the continuum limit, which have a maximal source-sink separation of 189 a .As an identity on the operator level, the PCAC relation is valid on every time sliceand, consequently, current quark masses exhibit a plateau, too. Contact terms lead todeviations from this plateau behaviour for small source-sink separations, as do cutoff effectsclose to the boundary. Following Ref. [23], we estimate the time slice extent of thesedeviations with the help of the same methods as for the effective meson masses above,i.e., modeling the corrections by exponential fits. For heavy-light and heavy-strange quarkmasses, the signal is lost in noise for large source-sink separations. Typical current quarkmass plateaux are displayed in the right part of Figure 1.In Appendix B we tabulate the ensemble-wise results for the meson masses and thedifferent bare current quark masses. In each table we list results before and after the massshifting procedure along the prescription in eq. (3.6). Accounting for these mass shiftsincreases the statistical uncertainties of the meson and bare quark masses considerably.This is due to the fact that the target for the mass shift, φ phys4 , is only known with about1 % precision, dominated by the uncertainty of the physical lattice scale. It should alsobe noted that the target value of φ phys4 is the same for all ensembles. For this reason theshifted masses from different gauge field ensembles cannot be considered uncorrelated. The PCAC definition (2.6) amounts to analysing the correlation functions f rs A P and f rs PP as described in the foregoing subsection, in order to arrive at numerical estimates for thecurrent quark masses referring to the various possible mass degenerate and non-degeneratecombinations of light, strange and heavy quarks. Assuming that a proper interpolationof these intermediate results involving the two different hopping parameter choices inthe charm region as a function of one of the two charmed meson masses (viz., m ¯D and m η c , advocated in Subsection 3.1 to fix the physical charm quark mass by their respectiveexperimental values) has been performed, one may define the renormalised charm quark12
20 40 60 80 100 120 x /a . . . . a m e ff h¯h, [43, 103]h¯s, [41, 101]h¯l, [40, 56] 0 20 40 60 80 100 120 x /a . . . . . . . a m r s h¯h , [32, 100]h¯s, [34, 110]h¯l, [41, 60]0 20 40 60 80 100 120 140 160 180 x /a . . . . . . a m e ff h¯h, [61, 158]h¯s, [65, 114]h¯l, [55, 94] 0 20 40 60 80 100 120 140 160 180 x /a . . . . . . . . a m r s h¯h , [53, 166]h¯s, [16, 132]h¯l, [16, 102] Figure 1:
Results from the extraction of meson and quark masses on the most chiral ensemble,D200 ( top ), and one of the ensembles at smallest lattice spacing, J501 ( bottom ). The labels denotethe quark content, where h corresponds to κ h of Table 2, and the intervals [ t min /a, t max /a ] indicatewhere the masses develop plateaux. Left : Effective meson masses that are employed in the tuning tothe physical charm quark mass.
Right : Bare current quark masses, employing improved derivativesin their lattice definition. Excessively fluctuating points based on heavy-light correlation functionsare not displayed. The uncertainties of the data points stemming from heavy-heavy correlationfunctions are too small to be visible on the chosen scale. mass from heavy-heavy correlation functions as m (c)R , c = m R , cc , (3.13)where c and c denote two distinct but mass degenerate quark flavours and the appropriaterenormalisation and improvement pattern is given in eq. (2.13). This definition alwaysprovides clean signals and therefore small statistical errors of the plateau averages. Thedependence of m (c)R , c on the light quark masses (encoded in its pion mass dependence) is13mall. Since the mass dependent cutoff effects are at the scale of the charm quark mass,we expect them to be sizeable, although effects of O( a ) are cancelled non-perturbatively inour lattice formulation of the theory.The size of these cutoff effects is reduced when heavy-light and heavy-strange correlationfunctions are employed from which estimators for the mass of the charm quark may thenbe obtained via m (l)R , c = 2 m R , lc − m R , ll , m (s)R , c = 2 m R , sc − m R , ss . (3.14)As evident in our data, the definitions in eq. (3.14) have a relevant dependence on thesea quark masses which, however, turns out to be significantly reduced in the associatedflavour-averaged charm quark mass m R , c = 23 m (l)R , c + 13 m (s)R , c , (3.15)because of the average sea quark mass being held constant along our trajectory in latticeparameter space towards the physical point. The signal of heavy-light and heavy-strangequark masses deteriorates for large source-sink separations. Accordingly, the statisticalerror on the charm quark mass definition m R , c is generally larger than on m (c)R , c .As a third option to arrive at a sensible prescription to compute the charm quark’smass, we adopt the ratio-difference method , eq. (2.19). Here, the first flavour is chosen tobe the heavy valence quark and the second flavour is one of the sea quarks. Adhering to thesame argument in favour of the linear combination in (3.15), we deduce the renormalisedcharm quark mass from the corresponding flavour average also in this case, i.e., m (rd)R , c = Z m r I , cl d I , cl r I , cl − r I , cs d I , cs r I , cs − ! , (3.16)with the improved ratio r I ,ij and difference d I ,ij already given in eqs. (2.17) and (2.18),respectively. The anticipated advantage of this estimator is that systematic uncertaintiesowing to an imprecise determination of the chiral point, which are naturally associatedwith subtracted quark masses, are cancelled in the quark mass from the ratio-differencemethod. As a central result of this paper, we wish to quote the renormalisation group invariant(RGI) value of the charm quark mass. For a given quark flavour, the RGI quark mass, M ,is defined through the formally exact expression M ≡ m R ( µ ) h b g ( µ ) i − d b exp ( − Z g ( µ )0 d g (cid:20) τ ( g ) β ( g ) − d b g (cid:21)) , (3.17)where it is assumed that renormalisation conditions are imposed at zero quark mass, whichsuggests calling those schemes mass independent. This implies ratios of renormalised quarkmasses for different flavours to become scale and scheme independent constants. In the sameway as the QCD Λ-parameter associated with the running of the renormalised coupling g R ( µ ), M belongs to the RGI quantities whose total dependence on the renormalisation14cale µ vanishes. In addition we recall that, even though perturbative expansions for theRG group functions β ( g ) and τ ( g ) exist, the Λ-parameter and the RGI quark masses aregenerically defined independent of perturbation theory and, particularly, the RGI quarkmasses M are renormalisation scheme independent. Therefore, we consider them idealfor comparisons with results from either experiments or other theoretical determinations.Conventions regarding the RG β - and τ -functions as well as their respective perturbativecoefficients b i and d i can be found, e.g., in Refs. [18, 61].The renormalised coupling g R ( µ ) and the continuum RG functions of the runningcoupling and quark mass entering eq. (3.17) were accurately obtained for the three-flavour theory at hand by applying the non-perturbative step scaling approach within theSchrödinger functional scheme [17, 62, 63]. Hence, we have all the ingredients available tomerge the above definition of M with one of the estimators for the renormalised charmquark mass introduced in the previous subsection and calculate the RGI charm quark massfrom it. Explicitly, on the basis of the expressions in eqs. (3.13), (3.15) and (3.16), theirrelation to the RGI mass then reads: M c = Mm R ( µ had ) m R , c ( µ had ) . (3.18)The renormalisation scale dependence of m R , c is inherited from the corresponding scaledependence of the respective Z -factor (see Section 2), which turns the underlying barequark mass into the renormalised one. Thanks to the determination of the universal,flavour independent ratio Mm R ( µ had ) = 0 . , (3.19)in the Schrödinger functional scheme for N f = 3 massless flavours at µ had = 233(8) MeV [17],the renormalisation factors and the ratio M/m R can be matched at the hadronic scale µ = µ had such that M c is readily returned. Further details concerning the non-perturbativecomputation of the running factor (3.19) can be found in Ref. [17] and references therein. After having extracted the renormalised quark masses as well as the relevant meson massesfrom every gauge field configuration ensemble, we proceed with a combined chiral andcontinuum limit extrapolation in order to obtain the value of the charm quark mass at thephysical point, defined by φ = φ phys2 and zero lattice spacing. All gauge field ensembleslisted in Table 1 are included in the analysis, except for H200 with a spatial extent of ∼ In the extrapolation, we include our different definitions (3.13),(3.15) and (3.16) of the renormalised charm quark mass, in conjunction with eq. (3.19)and made dimensionless by attaching the factor √ t . For the renormalisation constant Z , as well as the improvement coefficients b A − b P and b m entering these definitions, the Finite-volume effects for light current quark masses and the combined meson masses φ and φ onour set of ensembles were investigated in Refs. [23, 45] via a comparison of the results from the ensemblesH200 and N202. They were found to be significantly smaller than the available statistical accuracy. Weobserve the same to be true for the charmed observables studied here. In general, current quark massesand particularly masses of heavy mesons [64] are expected to be unaffected by moderate violations of theconditions that typically qualify lattice QCD ensembles as resembling an infinite volume situation. Z in Ref. [30], and found slightdifferences at finite lattice spacing but full agreement in the continuum limit. As for thesecond term in eq. (2.18) , we decided including it by recourse to the preliminary resultfor r m from Ref. [22]. However, our data reveals that neglecting this term only has aninsignificant effect on the final result.Eventually, after the joint extrapolation of √ t M c to the physical point, the finalvalue of the charm quark mass in physical units then simply stems from division by q t phys0 .For the charm quark mass in the continuum we assume the functional form √ t M continuumc ( φ , φ H ,
0) = c + c l φ + c c φ H , (4.1)which is guided by the following reasoning: Having corrected the chiral trajectory by meansof the mass derivatives in eq. (3.6), the light quark mass dependence is governed by the seapion mass only, whereas the kaon mass is implicitly fixed by the Tr[ M q ] = const . conditionrealised along the trajectory in lattice parameter space towards the physical point. Basedon the Gell-Mann-Oakes-Renner relation [65], we thus assume the leading contribution tobe proportional to φ . Moreover, we model the interpolation to the physical value of thebare charm quark mass via φ H = √ t m H , (4.2)where m H labels the mass of a charmed meson for which we either chose the flavouraveraged D meson mass, m ¯D , or the mass of the pseudoscalar charmonium meson, m η c (neglecting all quark-disconnected contributions). Motivated by heavy quark effectivetheory, we expect the charm quark mass to be linearly dependent on the meson mass ofchoice [54].For the continuum part we presume the leading cutoff effects to be of O( a ) as allrelevant renormalisation constants and improvement coefficients that propagate into ourcomputation are known non-perturbatively. In addition, we also allow for terms describingthe next two higher orders, O( a ) and O( a ). These cutoff effects may be quark massindependent or explicitly depend on the quark masses. Therefore, our general ansatz forthe lattice spacing dependence of the charm quark mass is parametrised by c M c ( φ , φ H , a ) = a t (cid:18) c + c φ (cid:19) + a (8 t ) / (cid:18) c + c φ (cid:19) + a (8 t ) (cid:18) c + c φ (cid:19) ; (4.3)all cutoff effects proportional to the light quark masses are already neglected here, becauseit turned out that these effects cannot be resolved in the fits. Furthermore, for now we alsoneglect terms proportional to mixed powers of a √ t and φ H , such as a (8 t ) / φ or a (8 t ) φ .We will demonstrate below that the inclusion of these terms (which a priori cannot beexcluded on theoretical grounds) has no effect on our final result. As argued at the end of Section 2, the third term in this equation can be safely ignored, because ¯ b m isonly of O( g ) in perturbation theory and poorly known non-perturbatively.
16o arrive at a combined model for the global fitting procedure, we can compose thechiral and continuum parts either linearly, by adding eqs. (4.1) and (4.3), viz. √ t M linearc ( φ , φ H , a ) = √ t M continuumc ( φ , φ H ,
0) + c M c ( φ , φ H , a ) , (4.4)or in a non-linear fashion by multiplication of the two: √ t M non − linearc ( φ , φ H , a ) = √ t M continuumc ( φ , φ H , × (1 + c M c ( φ , φ H , a )) ; (4.5)the c i in (4.1) and (4.3) are fit parameters. From now on we will refer to different modelparametrisations with the following labeling scheme:• The first entry is either ’l’ for the linear combination of the chiral and continuumparts or ’nl’ for the non-linear variant.• The second entry indicates the term which accounts for cutoff effects proportional to a . Here, ’1’ stands for the term parametrised by c , whereas ’2’ stands for the (heavymeson) mass dependent one, parametrised by c . When both pieces, proportional to c and c , are included, the entry is labelled by ’3’.• The third entry represents the contribution accounting for cutoff effects of O( a ).Here, ’0’ specifies fits without any contribution of this order, whereas ’1’ does so forfits with the c -term, ’2’ for fits with the c -term and ’3’ for fits with both of them.• The fourth entry labels the terms describing discretisation effects of O( a ). In analogyto the third entry, ’0’ labels no term of this order, ’1’ the term proportional to c , ’2’the one proportional to c and ’3’ both of them.• Finally, the fifth entry indicates whether a pion mass dependence is included in thefit, i.e., ’0’ marks models without and ’1’ models with the term parametrised by thecoefficient c l .As an example, the model labelled by (’l’, 3, 0, 1, 0) corresponds to the fit function √ t M ( l , , , , ( φ , φ H , a ) = c + c c φ H + a (cid:18) c t + c m (cid:19) + a (cid:18) c (8 t ) (cid:19) . (4.6)We always include at least one term describing cutoff effects of O( a ) and restrict themaximum number of terms related to the cutoff effects to three. In total this amounts to M = 108 different models, with which we attempt to fit our data.In order to estimate the optimal parameters for a given model (generically designated p i and f in the formula below), we employ an extended χ minimisation procedure, takinginto account the errors of both, the independent and the dependent variables [66], viz. χ = (cid:0) f ( p i ; ˜ x a ) − y a , ˜ x a − x a (cid:1) C − (cid:0) f ( p i ; ˜ x a ) − y a , ˜ x a − x a (cid:1) T . (4.7)The independent variables x a are promoted to fit parameters ˜ x a , which albeit stay con-strained to their initial values x a and thereby only indirectly affect the number of degreesof freedom. In our case, the full correlation matrix is ill-conditioned and eludes safe17nversion. For this reason we ignore all correlations and set the off-diagonal elements ofthe error covariance matrix C explicitly to zero.While minimising the uncorrelated χ yields reliable maximum likelihood estimatorsfor our parameters (see Ref. [67]), the minimal value of the uncorrelated χ cannot serveas a meaningful indicator for the goodness of fit, as we do not know the true number ofdegrees of freedom. In practice, we would thus obtain values of χ / ( N − k ) < χ ,the so-called ‘expected’ χ proposed and advocated for use in such situations in Ref. [60] ,which equals the expectation value of χ given normally distributed data with the fullcovariance matrix C . As explained in this reference, it qualifies to derive a correctedmeasure for the goodness of fit, χ = ( N − k ) χ χ , (4.8)in the sense that h χ i / ( N − k ) = 1 when the model describes the data. In order to address systematic effects in our computation of the charm quark massquantitatively, we classify our data in different categories w.r.t. the various choices andtechnical ingredients that characterise the analysis:• As estimator for the charm quark mass on the lattice, we pick the definitions m R , c , m (c)R , c or m (rd)R , c introduced in eqs. (3.13), (3.15) and (3.16).• For the bare current quark masses involved, either standard derivatives according toeq. (2.8) or its improved version in eq. (2.9) are used.• To fix the physical charm quark mass, we either employ the flavour-averaged D mesonmass, m ¯D , or the mass of the pseudoscalar charmonium groundstate, m η c (neglectingthe quark-disconnected contributions).• To ensure that cutoff effects are properly described by the fits, we either include allgauge field ensembles, i.e., covering a range of lattice spacings 0 . ≤ a fm ≤ . . ≤ a fm ≤ .
076 is imposed on the data, i.e., all ensembles at β = 3 . χ ,AIC = χ + 2 k , (4.9) Whereas Monte Carlo data from different gauge field ensembles is uncorrelated by construction, thefull covariance matrix does not assume a block diagonal form, but receives additional contributions inducedby common renormalisation and improvement parameters as well as the mass shifting procedure. Note that χ was already successfully applied before in a fit analysis of HQET correlators in thecontext of extracting semi-leptonic form factors [68]. .
00 0 .
01 0 .
02 0 .
03 0 . a / t . . . . . . √ t M c m R , c [ m η c , sd], (’l’, 1, 1, 2, 0) m R , c [ m η c , id], (’l’, 1, 2, 0, 0) m (c)R , c [ m η c , sd], (’nl’, 3, 2, 0, 0) m (c)R , c [ m η c , id], (’l’, 1, 2, 0, 0) .
75 3 .
80 3 .
85 3 .
90 3 .
95 4 .
00 4 .
05 4 .
10 4 . √ t m ¯D . . . . . . . . √ t M c m R , c [ m ¯D , sd], (’l’, 1, 1, 0, 0) m (c)R , c [ m ¯D , id], (’l’, 2, 2, 1, 0) m (rd)R , c [ m ¯D , sd], (’nl’, 2, 0, 2, 0) .
00 0 .
01 0 .
02 0 .
03 0 . a / t . . . . . . √ t M c m R , c [ m ¯D , sd], (’l’, 1, 1, 0, 0) m R , c [ m η c , sd], (’l’, 1, 1, 2, 0) m (rd)R , c [ m ¯D , id], (’nl’, 2, 0, 2, 0) m (rd)R , c [ m η c , id], (’l’, 1, 0, 1, 0) Figure 2:
Comparison of best fits according to the AIC in different categories. To map outtwo-dimensional representations of the fits, other dependencies than the ones illustrated in theplots have been projected to the physical point. The top panel reproduces the influence of utilisingstandard versus improved lattice derivatives to extract bare PCAC quark masses, while the middleone exhibits the clearly linear dependence of the charm quark mass estimators on the dimensionlessheavy meson mass, φ H = √ t m H , in case of m H = m ¯D . The bottom panel reflects the sensitivityof the joint global extrapolation results to the choice of φ H for fixing the physical charm point. Thelabeling convention as well as further details are given in the text. k is the number of parameters in the fit function.In Figure 2 we compare the best fits in different representative categories. The top partillustrates the effect of using improved derivatives in the PCAC masses on the definitionsof the charm quark mass. Starting with m R , c , we see that improvement of the derivativeleads to larger cutoff effects compared to the standard one. However, both variants agreevery well in the continuum limit. This is also the case for m (rd)R , c (not displayed in thefigure), but not so for m (c)R , c which is also shown in the top part of Figure 2. Here, thevariant with improved derivative coincides with the other results, while the extractionbased on the standard derivative tends to smaller values. Moreover, besides statisticaluncertainties being larger for the coarser lattice spacings, cutoff effects are expected tobe generically larger for doubly charmed PCAC masses. For these reasons we decide toexclude all determinations through m (c)R , c with the standard derivative in the final analysis.Regarding the lattice spacing dependence, we observe for all categories in our analysisthat terms of O( a ) or O( a ) have to be taken into account when considering ensembleswith a > .
06 fm. For finer resolutions, terms modelling a -contributions are sufficient todescribe the data reliably. The quadratic scaling regime therefore starts later than onecould have naively expected. Due to the purely non-perturbative inputs of renormalisationand improvement factors, lattice artefacts of O( a ) are absent.For all variants that we consider, fixing the physical value of the charm quark massworks completely satisfactory with a term linear in the (heavy) meson mass φ H = √ t m H .This is demonstrated in the middle part of Figure 2, which depicts the dependence ofthe three different charm quark mass definitions on m H = m ¯D . Employing m H = m η c instead to fix the physical charm quark mass, we find a very similar behaviour. Note thatfor β = 3 .
85 the two values chosen for κ c do not enclose the physical charm point; so inthis case we cannot interpolate, but have to rely on a slight extrapolation to the target.However, the perfectly linear behaviour reflected by the data suggests that this procedureis entirely legitimate.The bottom part of Figure 2 displays how our data and the respective joint chiral andcontinuum extrapolations depend on φ H = √ t m H that is chosen to fix the bare charmquark mass to its physical value. Compared to m H = m ¯D , the variants based on m η c ingeneral tend to slightly larger values at finite lattice spacing. However, this ambiguityvanishes in the continuum limit as expected.Since the fit distinguished as the best one by the AIC (or any other sensible measure)could in fact be also an outlier, we opt for performing a model averaging procedure withineach category, similar to what has been proposed in Ref. [71]. In analogy to this reference,we define the model average via h M c i = K X k =1 w k h M c i k , (4.10)with the respective model weights w k and k indexing the alltogether K = 108 models (i.e.,fit forms) within each of the 24 categories at disposal. For the systematic error arising A similar effect has been observed earlier in a computation of the D s meson decay constant in thequenched approximation [70]. σ M c = K X k =1 w k h M c i k − K X k =1 w k h M c i k ! . (4.11)For uniform weights, w k ≡ /K , this expression reduces to the variance of the differentmodel estimates. Yet instead of setting the weights uniformly, we rather set the w k foreach fit relying on the AIC, w AIC k = N exp (cid:18) −
12 AIC (cid:19) , (4.12)where the normalisation constant N is chosen to ensure P Kk =1 w AIC k = 1.In Figure 3 the outcome of the model averaging procedure is summarised for tworepresentative examples among the 24 categories. The individual circles correspond tothe result of fitting the data to the respective model, labelled along the horizontal axisaccording to scheme introduced above. Their opacity indicates the associated weight inthe average, the most opaque points having largest weight. The filled square and diamondmark the model averages with their statistical errors, while the shaded bands depict theirsystematic uncertainties emerging from the model averaging procedure. As a consequenceof the definitions of the renormalised charm quark mass, eqs. (3.13), (3.15) and (3.16),which were suitably adapted to the position of the CLS ensembles in physical parameterspace and the approach towards the chiral limit ensuing from it, we expect the (sea) pionmass dependence of M c to be small. In fact, the best fits in all 24 categories (and also themajority of fits with a relevant weight) do not incorporate a pion mass dependence.As can be inferred from the top part of Figure 3, almost all fits with non-negligibleweight agree very well, which is reflected by the small systematic error. In this case, singlingout only the best fit would give nearly an identical result compared to the model average.In the bottom part of Figure 3, the best fit is (’nl’, 2, 1, 1, 0) and differs from the finalmodel average by ∼ σ . However, owing to the very construction of this average, thelarge number of models that describe the data reasonably, outweigh the few outliers. Thesomewhat wider spread of the fits with more significant weight is revealed by the largersystematic model averaging error in this case.Figure 4 displays the mean values and errors from the model averages within each ofthe 24 categories specified above. In principle, the result in every single category representsa valid determination of the charm quark mass, which means that for infinite statisticswe would expect all 24 determinations to coincide perfectly. Hence, to decrease the finalerror, one could even consider combined fits of different categories. However, as we willdiscuss in the next section, our error is not dominated by the statistical error arising fromthe underlying Monte Carlo data, but by the propagated error from external inputs. Forthis reason we have decided to conservatively take the spread of the individual categories’results as additional measure for the systematic error of our calculation.For our final result of M c we chose the weighted average over all categories in Figure 4,except (as argued before) for the ones with m (c)R , c and the standard derivative. The weightsare given by the inverse of the respective errors (i.e., statistical and systematic ones addedin quadrature). From the weighted standard deviation (4.11) of the model averages withinthese categories we estimate the size of the systematic error in our determination of the21 o d e l a v e r ag e ( ’l’, , , , )( ’l’, , , , )( ’l’, , , , )( ’l’, , , , )( ’l’, , , , )( ’l’, , , , )( ’l’, , , , )( ’l’, , , , )( ’l’, , , , )( ’l’, , , , )( ’l’, , , , )( ’l’, , , , )( ’l’, , , , )( ’l’, , , , )( ’l’, , , , )( ’l’, , , , )( ’l’, , , , )( ’l’, , , , )( ’l’, , , , )( ’l’, , , , )( ’l’, , , , )( ’l’, , , , )( ’l’, , , , )( ’l’, , , , )( ’l’, , , , )( ’l’, , , , )( ’l’, , , , )( ’ n l’, , , , )( ’ n l’, , , , )( ’ n l’, , , , )( ’ n l’, , , , )( ’ n l’, , , , )( ’ n l’, , , , )( ’ n l’, , , , )( ’ n l’, , , , )( ’ n l’, , , , )( ’ n l’, , , , )( ’ n l’, , , , )( ’ n l’, , , , )( ’ n l’, , , , )( ’ n l’, , , , )( ’ n l’, , , , )( ’ n l’, , , , )( ’ n l’, , , , )( ’ n l’, , , , )( ’ n l’, , , , )( ’ n l’, , , , )( ’ n l’, , , , )( ’ n l’, , , , )( ’ n l’, , , , )( ’ n l’, , , , )( ’ n l’, , , , )( ’ n l’, , , , )( ’ n l’, , , , ) M c i n M e V M c i n M e V Figure 3:
Model averaging procedure for two representative categories: m R , c with m H = m ¯D ,standard lattice derivative and all lattice spacings included (top panel); and m (rd)R , c with m H = m η c ,improved derivative and all lattice spacings included (bottom panel). Each circle corresponds to theresult of the model indicated along the horizontal axis. The opacity of the circles is proportionalto their weights contributing to the average. The filled square and diamond give the respectivemodel averages with statistical error according to eq. (4.10), while the transparent bands indicatethe systematic errors arising from the modeling procedure according to eq. (4.11), both of whichare added in quadrature. Fits with pion mass dependence are not shown, as their contributionto the model averages is subdominant. Note that the model averages obtained as visualised herecorrespond to one data point each in the subsequent Figure 4. M c i n M e V Final estimate m H = m ¯D m H = m η c Standard derivativeImproved derivative0 . ≤ a fm ≤ . . ≤ a fm ≤ . m R , c m (c)R , c m (rd)R , c Figure 4:
Model averages in each category as evaluated through the AIC weighting prescription.Results for the m (c)R , c variant with the standard derivative, depicted as transparent points, areexcluded from our final estimates as explained in the text. The shaded region indicates the estimateof the systematic error obtained by the weighted standard deviation of the model averages over allcategories. (The legend recalls the key characteristics of every category.) charm quark mass. However, we have convinced ourselves that an unweighted averageleads to a number which is indistinguishable from the weighted average in relation to thetotal error for the final result quoted below.In order to check the robustness of our analysis method, we have also induced differentweights in the model averaging procedure, devised to set different penalties on the numberof model parameters and thereby to favour over- or underfitting opposed to the standardprocedure. Specifically, the Bayesian information criterion [72] as well as the correctedAIC [73] were tried, which give higher penalties for models with more fit parameters,while basing the weights on just χ introduces no penalty at all. The resultingvariation on our final estimate of M c when using these different weights is by far smallerthan our systematic error as assessed above. Furthermore, as already mentioned earlier,we investigated the effect of allowing for terms such as a √ t ( am H ) or a t ( am H ) in thegeneric fit ansatz for the joint chiral and continuum limit extrapolation, thus increasingthe number of possible models in each category to K = 168. Again, also this extension ofthe number of models has a completely negligible effect on the final result. In summary, wetherefore conclude that grounding our analysis on a large number of models and categoriesvery consistently yields a reliable result that is fairly independent of human bias, becausethe manual input is only limited to a theoretically well justified pre-selection of consideredcategories and models. 23 able 3: Error budget of our computation, quantified by the contributions to the squared errorof our final estimate for M c . We group them into the three main categories ‘Statistical error’,‘Renormalisation group running’ and ‘Systematic uncertainty’, which comprise the final (total)error, where the listed individual subcategories indicate how these uncertainties are distributedamong the main categories. Statistical error 44.8% correlation functions 2.2%O( a ) improvement 14.2%renormalisation 8.7%scale setting 19.7%physical meson masses < 0.1% Renormalisation group running 45.9%Systematic error 9.3% model selection 1.1%input data variation 8.2%
Following the data analysis presented in the foregoing sections, we quote for the RGI charmquark mass from our lattice QCD computation as final result M c ( N f = 3) = 1486(14)(14)(6) MeV = 1486(21) MeV , (5.1)which constitutes an unambiguous result in (the continuum limit of) the N f = 2 + 1 theorystudied here. The first error is statistical, the second stems from the factor convertingrenormalised into RGI quark masses, available from Ref. [17], while the third one representsour estimate of systematic effects.In Table 3 we summarise the error budget of our computation. The most dominantcontribution to the total error originates from the (universal) conversion factor M/m R ( µ had )that connects renormalised quark masses at finite scale via non-perturbative renormalisationgroup running to their RGI counterparts (cf. Subsection 3.4) and was obtained for thethree-flavour theory in Ref. [17]. The next relevant contributions come from the scalesetting of Ref. [23], and the O( a ) improvement b -coefficients multiplying the quark massdependent pieces, non-perturbatively obtained in Ref. [21]. For the systematic uncertainty,which is made up of two contributions, the one from the spread of outcomes of physicalpoint extrapolations within different categories (i.e., characterised by variations of theinput data as outlined in Subsection 4.1) outweighs the other one from the model averagingprocedure itself (i.e., involving the various fit forms). Note that the statistical error ofthe correlation functions entering the individual definitions of the charm quark mass is byfar subleading in our calculation and that the uncertainties on the experimental values ofthe meson masses, which feed into the analysis to set the physical point, are completelyirrelevant for the total error of our central result (5.1).As the final uncertainty of our charm quark mass determination is completely domi-nated by external input, there is hardly any room for a significant gain in precision, even ifthe underlying data base of gauge field configuration ensembles on the Tr[ M q ] = const . trajectory in lattice parameter space were moderately larger. One key prerequisite for24n improved overall precision on the RGI (charm) quark mass along the methodologypresented in this work would be a more accurate knowledge of the non-perturbative quarkmass running factor. A promising step in that direction is the idea of a renormalisationstrategy proposed in Ref. [74] that exploits the decoupling of heavy quarks from low-energyphysics. Another substantial refinement in precision may be anticipated from an update ofthe physical scale setting parameter, q t phys0 , as it is currently being pursued within theCLS effort (see Refs. [75, 76]).For comparison to other results (from lattice as well as non-lattice QCD determinations)as well as later reference, it is customary to convert our result to the four-flavour theorythrough the well-established prescription relying on perturbative decoupling in the MSscheme. This procedure is sketched in some detail in Appendix A, where we also providevalues of the charm quark mass in the MS scheme at conventional energy scales. Particularly,for the four-flavour RGI mass we arrive at: M c ( N f = 4 , Λ MeV = 1548(23) MeV . (5.2)Its error splits into the part derived from the uncertainty of the three-flavour RGI charmquark mass above, to which the subdominant part invoked by Λ (3)MS = 341(12) MeV [77](that we employ for the QCD Λ-parameter entering the perturbative renormalisation grouprunning) is added in quadrature. We have calculated the charm quark mass in QCD with N f = 2 + 1 dynamical quarkflavours, on the basis of a large set of gauge field configuration ensembles generatedby the Coordinated Lattice Simulations (CLS) initiative [4, 5] and employing a partiallyquenched setup, consisting of non-perturbatively O( a ) improved Wilson quarks in the seasupplemented by a quenched charm quark in the valence sector. Bare current quark masseswere extracted from the PCAC relation. Notable features of our determination are: five(very) fine lattice spacings from 0 .
087 fm down to less than 0 .
04 fm, fully non-perturbativeO( a ) improvement of the quark mass and the use of the non-perturbative renormalisationgroup running, available from Ref. [17], in order to arrive at the RGI (renormalisationgroup invariant, i.e., scheme and scale independent) estimate of the charm quark massas central result. Moreover, our analysis exploits three different definitions of the quarkmass on the lattice, which also allows systematic effects to be investigated in a controlledway. To ensure reliability and robustness of our result that, for given possible variations inthe analysis, includes a large enough (albeit theoretically well-founded) number of modelfunctions for the joint chiral and continuum extrapolation whose individual contributionsare correctly weighted by their fit to the actual data, we have adopted a model averagingprocedure inspired by the prescription proposed in Ref. [71]. In this regard, our finalresult can therefore be considered as model-independent and its error as accounting for allsystematic effects of the present computation in the studied N f = 2 + 1 approximation ofQCD.The final numerical outcome of our work is to be taken from eqs. (5.1) and (5.2) ofthe previous section, where we quote the RGI charm quark mass in the three-flavour (i.e., N f = 2 + 1) and four-flavour theories, respectively, the latter being obtained by virtue of25erturbative decoupling relations. It is worth emphasising that the precision of our resultcan only be substantially improved by reducing the uncertainties on external quantities(most notably, the renormalised-to-RGI quark mass conversion factor and the lattice scalein physical units) which enter our calculation. Furthermore, while the joint extrapolationsof the charm quark mass to the physical point exhibit, in line with O( a ) improvement, nocorrections linearly, but quadratically and of higher order in the lattice spacing, they arenot able to resolve any dependence on the sea pion mass; this stands in distinct contrastto the calculation of light quark masses on a subset of the gauge configuration ensemblesat hand, reported in Ref. [45], for which a decent decrease of errors can be expected fromaccommodating more chiral ensembles (especially with fine lattice spacings) and morestatistics at light quark masses.As soon as the relative error is smaller than about 1%, systematic uncertainties dueto neglecting the charm quark in the sea as well as isospin splitting and QED effects maybecome noticeable. However, as apparent from Figure 5 where we compare our result withother lattice QCD ones from the literature, no significant difference between 2 + 1 and2 + 1 + 1 flavour results can be unravelled at the current level of precision. For instance,the impact of the inclusion of quenched QED was found to be around 0 .
2% in Ref. [78],which is well below our final uncertainty.After all, in Appendix A, we explain how we translated eqs. (5.1) and (5.2) to associatedcharm quark mass values in the MS scheme, present our results at customary energy scales,and confront them with those of other lattice QCD collaborations. Here, we compare ourfour-flavour RGI mass estimate M c ( N f = 4) = 1548(23) MeV to the current FLAG averagefor calculations based on N f = 2 + 1 dynamical quarks [15], which is composed of thenumbers from Refs. [7–9] and reads: M FLAGc ( N f = 4) = 1529(17) MeV . (6.1)Within their 1 σ -errors, both are consistent with each other. Acknowledgements.
We thank Mattia Bruno, Tomasz Korzec, Stefan Schaefer and RainerSommer, as well as Gunnar Bali, Sjoerd Bouma, Sara Collins and Wolfgang Söldner, for helpfuldiscussions and the RQCD Collaboration for the fruitful collaboration in our joint project oncharmed decay constants. We are also indebted to Patrick Fritzsch, Carlos Pena and AnastassiosVladikas for a critical reading of an earlier version of the manuscript and their valuable comments.We are grateful to our colleagues in the CLS initiative for producing the gauge field configurationensembles used in this study. This work is supported by the Deutsche Forschungsgemeinschaft(DFG) through the Research Training Group
GRK 2149: Strong and Weak Interactions – fromHadrons to Dark Matter . We acknowledge the computer resources provided by the WWU IT,formerly ‘Zentrum für Informationsverarbeitung (ZIV)’, of the University of Münster (PALMA-IIHPC cluster) and thank its staff for support. Perturbative quark mass running to conventional scales
In view of further applications or reference in phenomenological studies, the charm quarkmass is often quoted in the MS renormalisation scheme and the four-flavour theory. In thecase of heavy quarks, in particular, it is also common practice to do this in form of thescale invariant mass, which amounts to implicitly fixing the renormalisation (i.e., energy)scale of the running MS quark mass at its own value.For the purpose of this conversion, starting from the RGI mass value, we first usethe result Λ (3)MS = 341(12) MeV for the QCD Λ-parameter computed by the ALPHACollaboration in Ref. [77], together with perturbative running in the MS scheme and thethree-flavour theory, to obtain the scale invariant charm quark mass, m MSc ( µ = m MSc ; N f =3); for the numerical evaluation, we proceed along the lines of Appendix A of Ref. [79].Next, at the energy scale defined by the scale invariant charm quark mass, one has toemploy (perturbative) decoupling relations, which account for the effect of the changein the number of active quark flavours when passing across the respective quark massthresholds while running with energy, in order to consistently arrive at a result in thefour-flavour theory; for this step, we utilise the RunDec software package described inRefs. [80–82]. Then one can perturbatively run this number in the four-flavour theory,either to the scale invariant quark mass (which gets slightly shifted by the decouplingprocess) or upwards to the conventional energy scale µ = 3 GeV. In this way we provideresults, for both 4-loop and 5-loop perturbative accuracy of the β -function and the quarkmass anomalous dimension , and find: m MSc ( µ = m MSc ; N f = 4 , Λ MeV = 1292(19) MeV , (A.1) m MSc ( µ = m MSc ; N f = 4 , Λ MeV = 1296(19) MeV ; (A.2)and m MSc ( µ = 3 GeV; N f = 4 , Λ MeV = 1005(16) MeV , (A.3) m MSc ( µ = 3 GeV; N f = 4 , Λ MeV = 1007(16) MeV , (A.4)where the first error is the one from our determination discussed in the main text, whereasthe second error stems from the uncertainty of Λ (3)MS .After the switch to the four-flavour theory, we may employ perturbative running againin order to eventually obtain the four-flavour RGI charm quark mass, for which we quote M c ( N f = 4 , Λ MeV = 1548(23) MeV , (A.5)and the 5-loop value already stated in eq. (5.2). Note that the error contribution from Λ (3)MS becomes negligible in the four-flavour RGI limit and that the 4-loop and 5-loop resultsagree up to less than 1 MeV. We thus conclude that any ambiguities, which arise from theperturbative running and decoupling in the MS scheme, are much smaller than the totalerror of our final result.In Figure 5 we compare our result for m MSc ( m MSc ) in the four-flavour theory with thosefrom other lattice QCD determinations based on gauge field configurations with 2 + 1 and Thanks to this high-order accuracy, perturbation theory that has to be invoked within any conversionto the conventionally cited MS scheme looks very well behaved. .
24 1 .
26 1 .
28 1 .
30 1 .
32 1 .
34 1 .
36 1 . m MSc ( m MSc ) in GeV FLAG 2+1+1 averageHPQCD 18FNAL/MILC/TUMQCD 18HPQCD 14AETM 14AETM 14FLAG 2+1 averageJLQCD 16 χ QCD 14HPQCD 10PDGthis work N f = 2 + 1 + 1 N f = 2 + 1 Figure 5:
Comparison of our result (orange point) for m MSc ( m MSc ) in the four-flavour theory with1 . N f = 2+1+1 and N f = 2+1 FLAG averages [15](grey vertical bands), as well as the lattice results that entered these averages: HPQCD 18 [14],FNAL/MILC/TUMQCD 18 [13], HPQCD 14A [12], ETM 14A [11], ETM 14 [10], JLQCD 16 [9], χ QCD 14 [8], HPQCD 10 [7]. The upper part of the figure displays lattice results originating fromsimulations with four quark flavours in the sea (marked by diamonds), whereas the lower part refersto those obtained on gauge configurations with three dynamical quark flavours (marked by squares). σ uncertainties. Giventhat the underlying strategies to extract the physical quark mass, as well as the discretisedactions, renormalisation procedures and further potential systematics, differ significantlyamong these calculations, this overall agreement appears very reassuring.Let us close this appendix with a general remark on the renormalisation strategywhich we follow in this work. Indeed, one might argue that it is not ideally suited whenone is simply interested in the charm quark mass in the MS scheme at a scale of a fewGeV: namely, as a consequence of first running non-perturbatively (in the intermediateSchrödinger functional scheme) up to an energy scale of about M W deeply in the per-turbative regime, from which the scale evolution of the renormalised mass can safelybe continued perturbatively to infinite energy where the RGI mass is defined, and then,again perturbatively, down to a few GeV in the MS scheme, the errors from both runningprocedures dominate our uncertainty on m MSc ( µ = m MSc ) by far. Based on this observation,one could come up with the idea to only run non-perturbatively to the scale of choice, e.g.,2 GeV, and to perform the perturbative conversion to the MS scheme directly at this scale.However, apart from the systematic uncertainties that are connected with the perturbativeconversion at low-energy scales, the statistical error would not be decreased noticeably,since the dominant contribution to the uncertainty of the relevant factor encoding thenon-perturbative running of the renormalised quark mass (analogous to (3.19)) is alreadyaccumulated in the low-energy regime, i.e., below 2 GeV. Apart from these conceptual28isadvantages of the MS mass within our approach, we also prefer the RGI quark massbecause of its a priori genuinely non-perturbative nature as well as its renormalisationscheme and scale independence and, therefore, distinctly more controllable systematicerror. 29
Tables
Table 4:
Overview of the results for the gluonic quantity t /a as well as the combinations φ and φ defined in eqs. (3.1) and (3.3) from pion and kaon masses, together with the associated hoppingparameters. For every gauge field configuration ensemble we list the numbers before the quarkmass shift (see Section 3) in the rows with labels ‘id’ and the corresponding ones after the shift to φ phys4 in the rows without labels.id t /a φ φ κ l κ s H101 2 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . able 5: Overview of the results for the flavour-averaged heavy-light and heavy-heavy pseudoscalarmeson masses, referring to the two hopping parameter choices in the charm region given in Table 2.As before, for every ensemble we list the numbers before the shift in the rows with labels and thecorresponding ones after the shift to φ phys4 in the rows without labels.id am ¯ D am ¯ D am η c , am η c , H101 0 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . able 6: Overview of the results for light-light, strange-strange and heavy-light PCAC quarkmasses from standard lattice derivatives, where the heavy quark hopping parameter choices in thecharm region refer to Table 2. As before, for every ensemble we list the numbers before the shift inthe rows with labels and the corresponding ones after the shift to φ phys4 in the rows without labels.Analogous results from improved lattice derivatives are collected in Table 8.id am ll am ss am h l am h l H101 0 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . able 7: Overview of the results for heavy-strange and heavy-heavy PCAC quark masses fromstandard lattice derivatives, referring to the heavy quark hopping parameter choices given in Table 2.As before, for every ensemble we list the numbers before the shift in the rows with labels and thecorresponding ones after the shift to φ phys4 in the rows without labels. Analogous results fromimproved lattice derivatives are collected in Table 9.id am h s am h s am h h am h h H101 0 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . able 8: Overview of the results for light-light, strange-strange and heavy-light PCAC quarkmasses from improved lattice derivatives, where the heavy quark hopping parameter choices in thecharm region refer to Table 2. As before, for every ensemble we list the numbers before the shift inthe rows with labels and the corresponding ones after the shift to φ phys4 in the rows without labels.Analogous results from standard lattice derivatives are collected in Table 6.id am ll am ss am h l am h l H101 0 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . able 9: Overview of the results for heavy-strange and heavy-heavy PCAC quark masses fromimproved lattice derivatives, referring to the heavy quark hopping parameter choices given in Table2. As before, for every ensemble we list the numbers before the shift in the rows with labels andthe corresponding ones after the shift to φ phys4 in the rows without labels. Analogous results fromstandard lattice derivatives are collected in Table 7.id am h s am h s am h h am h h H101 0 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . eferences [1] M. Antonelli et al., Flavor Physics in the Quark Sector , Phys. Rept. (2010) 197,[ ].[2] J. R. Andersen et al.,
Handbook of LHC Higgs Cross Sections: 3. Higgs Properties , .[3] A. A. Petrov, S. Pokorski, J. D. Wells and Z. Zhang, Role of low-energy observables inprecision Higgs boson analyses , Phys. Rev. D (2015) 073001, [ ].[4] M. Bruno et al., Simulation of QCD with N f = + , JHEP (2015) 043, [ ].[5] D. Mohler, S. Schaefer and J. Simeth, CLS 2+1 flavor simulations at physical light- andstrange-quark masses , EPJ Web Conf. (2018) 02010, [ ].[6] S. Cali, F. Knechtli and T. Korzec,
How much do charm sea quarks affect the charmoniumspectrum? , Eur. Phys. J. C (2019) 607, [ ].[7] C. McNeile, C. Davies, E. Follana, K. Hornbostel and G. Lepage, High-Precision c and bMasses, and QCD Coupling from Current-Current Correlators in Lattice and ContinuumQCD , Phys. Rev. D (2010) 034512, [ ].[8] Y.-B. Yang et al., Charm and strange quark masses and f D s from overlap fermions , Phys.Rev. D (2015) 034517, [ ].[9] K. Nakayama, B. Fahy and S. Hashimoto, Short-distance charmonium correlator on the latticewith Möbius domain-wall fermion and a determination of charm quark mass , Phys. Rev. D (2016) 054507, [ ].[10] N. Carrasco et al., Up, down, strange and charm quark masses with N f = 2+1+1 twistedmass lattice QCD , Nucl. Phys. B (2014) 19, [ ].[11] C. Alexandrou, V. Drach, K. Jansen, C. Kallidonis and G. Koutsou,
Baryon spectrum with N f = 2 + 1 + 1 twisted mass fermions , Phys. Rev. D (2014) 074501, [ ].[12] B. Chakraborty, C. Davies, B. Galloway, P. Knecht, J. Koponen, G. C. Donald et al., High-precision quark masses and QCD coupling from n f = 4 lattice QCD , Phys. Rev. D (2015) 054508, [ ].[13] A. Bazavov et al., Up-, down-, strange-, charm-, and bottom-quark masses from four-flavorlattice QCD , Phys. Rev. D (2018) 054517, [ ].[14] A. Lytle, C. Davies, D. Hatton, G. Lepage and C. Sturm, Determination of quark masses from n f = 4 lattice QCD and the RI-SMOM intermediate scheme , Phys. Rev. D (2018) 014513,[ ].[15] (Flavour Lattice Averaging Group) , S. Aoki et al., FLAG Review 2019: Flavour LatticeAveraging Group (FLAG) , Eur. Phys. J. C (2020) 113, [ ].[16] S. Weinberg, New approach to the renormalization group , Phys. Rev. D (1973) 3497.[17] I. Campos, P. Fritzsch, C. Pena, D. Preti, A. Ramos and A. Vladikas, Non-perturbative quarkmass renormalisation and running in N f = 3 QCD , Eur. Phys. J. C (2018) 387,[ ].[18] S. Capitani, M. Lüscher, R. Sommer and H. Wittig, Non-perturbative quark massrenormalization in quenched lattice QCD , Nucl. Phys. B (1999) 669, [ hep-lat/9810063 ].[Erratum: Nucl. Phys. B (2000) 762].
19] J. Bulava, M. Della Morte, J. Heitger and C. Wittemeier,
Non-perturbative improvement ofthe axial current in N f =3 lattice QCD with Wilson fermions and tree-level improved gaugeaction , Nucl. Phys. B (2015) 555, [ ].[20] M. Dalla Brida, T. Korzec, S. Sint and P. Vilaseca,
High precision renormalization of theflavour non-singlet Noether currents in lattice QCD with Wilson quarks , Eur. Phys. J. C (2019) 23, [ ].[21] G. M. de Divitiis, P. Fritzsch, J. Heitger, C. C. Köster, S. Kuberski and A. Vladikas, Non-perturbative determination of improvement coefficients b m and b A − b P and normalisationfactor Z m Z P /Z A with N f = 3 Wilson fermions , Eur. Phys. J. C (2019) 797, [ ].[22] J. Heitger, F. Joswig, P. L. J. Petrak and A. Vladikas, The ratio of the renormalisationparameters of the flavour singlet and non-singlet scalar densities in N f = 3 lattice QCD withWilson quarks , in preparation.[23] M. Bruno, T. Korzec and S. Schaefer, Setting the scale for the CLS flavor ensembles , Phys. Rev. D (2017) 074504, [ ].[24] J. Heitger, F. Joswig and S. Kuberski, Towards the determination of the charm quark mass on N f = 2 + 1 CLS ensembles , PoS
LATTICE2019 (2019) 092, [ ].[25] M. Bruno, I. Campos, J. Koponen, C. Pena, D. Preti, A. Ramos et al.,
Light and strangequark masses from N f = 2 + 1 simulations with Wilson fermions , PoS
LATTICE2018 (2019) 220, [ ].[26] T. Bhattacharya, R. Gupta, W. Lee, S. R. Sharpe and J. M. S. Wu,
Improved bilinears inlattice QCD with non-degenerate quarks , Phys. Rev. D (2006) 034504, [ hep-lat/0511014 ].[27] G. de Divitiis and R. Petronzio, Nonperturbative renormalization constants on the lattice fromflavor nonsinglet Ward identities , Phys. Lett. B (1998) 311, [ hep-lat/9710071 ].[28] M. Guagnelli, R. Petronzio, J. Rolf, S. Sint, R. Sommer and U. Wolff,
Nonperturbative resultsfor the coefficients b m and b A − b P in O( a ) improved lattice QCD , Nucl. Phys. B (2001)44, [ hep-lat/0009021 ].[29] S. Dürr, Z. Fodor, C. Hoelbling, S. D. Katz, S. Krieg, T. Kurth et al.,
Lattice QCD at thephysical point: Simulation and analysis details , JHEP (2011) 148, [ ].[30] J. Heitger, F. Joswig and A. Vladikas, Ward identity determination of Z S /Z P for N f = 3 lattice QCD in a Schrödinger functional setup , Eur. Phys. J. C (2020) 765, [ ].[31] J. Heitger and J. Wennekers, Effective heavy light meson energies in small volume quenchedQCD , JHEP (2004) 064, [ hep-lat/0312016 ].[32] P. Fritzsch, J. Heitger and N. Tantalo, Non-perturbative improvement of quark massrenormalization in two-flavour lattice QCD , JHEP (2010) 074, [ ].[33] P. Korcyl and G. S. Bali, Non-perturbative determination of improvement coefficients usingcoordinate space correlators in N f = 2 + 1 lattice QCD , PoS
LATTICE2016 (2016) 190,[ ].[34] P. Korcyl and G. S. Bali,
Non-perturbative determination of improvement coefficients usingcoordinate space correlators in N f = 2 + 1 lattice QCD , Phys. Rev. D (2017) 014505,[ ].[35] G. S. Bali, E. E. Scholz, J. Simeth and W. Söldner, Lattice simulations with N f = 2 + 1 improved Wilson fermions at a fixed strange quark mass , Phys. Rev. D (2016) 074501,[ ].
36] M. Lüscher and P. Weisz,
On-Shell Improved Lattice Gauge Theories , Commun. Math. Phys. (1985) 59. [Erratum: Commun. Math. Phys. (1985) 433].[37] B. Sheikholeslami and R. Wohlert, Improved Continuum Limit Lattice Action for QCD withWilson Fermions , Nucl. Phys. B (1985) 572.[38] M. Lüscher and S. Schaefer,
Lattice QCD without topology barriers , JHEP (2011) 036,[ ].[39] M. Lüscher and S. Schaefer, Lattice QCD with open boundary conditions and twisted-massreweighting , Comput. Phys. Commun. (2013) 519, [ ].[40] M. Lüscher and F. Palombi,
Fluctuations and reweighting of the quark determinant on largelattices , PoS
LATTICE2008 (2008) 049, [ ].[41] A. D. Kennedy, I. Horvath and S. Sint,
A New Exact Method for Dynamical FermionComputations with Non-Local Actions , Nucl. Phys. Proc. Suppl. (1999) 834,[ hep-lat/9809092 ].[42] M. A. Clark and A. D. Kennedy, Accelerating dynamical fermion computations using therational hybrid Monte Carlo (RHMC) algorithm with multiple pseudofermion fields , Phys. Rev.Lett. (2007) 051601, [ hep-lat/0608015 ].[43] M. Bruno, J. Finkenrath, F. Knechtli, B. Leder and R. Sommer, Effects of Heavy Sea Quarksat Low Energies , Phys. Rev. Lett. (2015) 102001, [ ].[44] S. Cali, F. Knechtli and T. Korzec,
Comparison between models of QCD with and withoutdynamical charm quarks , PoS
LATTICE2018 (2018) 082, [ ].[45] M. Bruno, I. Campos, P. Fritzsch, J. Koponen, C. Pena, D. Preti et al.,
Light quark masses in N f = 2 + 1 lattice QCD with Wilson fermions , Eur. Phys. J. C (2020) 169, [ ].[46] M. Lüscher, Properties and uses of the Wilson flow in lattice QCD , JHEP (2010) 071,[ ]. [Erratum: JHEP (2014) 092].[47] (Flavour Lattice Averaging Group) , S. Aoki et al., Review of lattice results concerninglow-energy particle physics , Eur. Phys. J. C (2017) 112, [ ].[48] T. Korzec, https://github.com/to-ko/mesons .[49] M. Bruno, P. Korcyl, T. Korzec, S. Lottini and S. Schaefer, On the extraction of spectralquantities with open boundary conditions , PoS
LATTICE2014 (2014) 089, [ ].[50] G. M. de Divitiis, R. Petronzio and N. Tantalo,
Distance preconditioning for lattice Diracoperators , Phys. Lett. B (2010) 157, [ ].[51] S. Collins, K. Eckert, J. Heitger, S. Hofmann and W. Söldner,
Charmed pseudoscalar decayconstants on three-flavour CLS ensembles with open boundaries , PoS
LATTICE2016 (2017)368, [ ].[52] P. de Forcrand, M. Garcia Perez, H. Matsufuru, A. Nakamura, I. Pushkina, I.-O. Stamatescuet al.,
Contribution of disconnected diagrams to the hyperfine splitting of charmonium , JHEP (2004) 004, [ hep-lat/0404016 ].[53] A. F. Falk and M. Neubert, Second order power corrections in the heavy quark effective theory.1. Formalism and meson form-factors , Phys. Rev. D (1993) 2965, [ hep-ph/9209268 ].[54] M. Neubert, Heavy quark symmetry , Phys. Rept. (1994) 259, [ hep-ph/9306320 ].[55] S. Collins and W. Söldner, private communication.[56] S. Collins, K. Eckert, J. Heitger, S. Hofmann and W. Söldner,
Leptonic decay constants forD-mesons from 3-flavour CLS ensembles , EPJ Web Conf. (2018) 13019, [ ].
57] U. Wolff,
Monte Carlo errors with less errors , Comput. Phys. Commun. (2004) 143,[ hep-lat/0306017 ]. [Erratum: Comput. Phys. Commun. (2007) 383].[58] A. Ramos,
Automatic differentiation for error analysis of Monte Carlo data , Comput. Phys.Commun. (2019) 19, [ ].[59] S. Schaefer, R. Sommer and F. Virotta,
Critical slowing down and error analysis in latticeQCD simulations , Nucl. Phys. B (2011) 93, [ ].[60] M. Bruno and R. Sommer, in preparation.[61] M. Della Morte, R. Hoffmann, F. Knechtli, J. Rolf, R. Sommer, I. Wetzorke et al.,
Non-perturbative quark mass renormalization in two-flavor QCD , Nucl. Phys. B (2005)117, [ hep-lat/0507035 ].[62] M. Dalla Brida, P. Fritzsch, T. Korzec, A. Ramos, S. Sint and R. Sommer,
Determination ofthe QCD Λ -parameter and the accuracy of perturbation theory at high energies , Phys. Rev.Lett. (2016) 182001, [ ].[63] M. Dalla Brida, P. Fritzsch, T. Korzec, A. Ramos, S. Sint and R. Sommer,
Slow running ofthe Gradient Flow coupling from 200 MeV to 4 GeV in N f = 3 QCD , Phys. Rev. D (2017)014507, [ ].[64] G. Colangelo, A. Fuhrer and S. Lanz, Finite volume effects for nucleon and heavy mesonmasses , Phys. Rev. D (2010) 034506, [ ].[65] M. Gell-Mann, R. J. Oakes and B. Renner, Behavior of current divergences under SU × SU , Phys. Rev. (1968) 2195.[66] P. T. Boggs and J. E. Rogers,
Orthogonal distance regression , tech. rep., National Institute ofStandards and Technology, Gaithersburg, MD, 1989. 10.6028/NIST.IR.89-4197.[67] C. Michael,
Fitting correlated data , Phys. Rev. D (1994) 2616, [ hep-lat/9310026 ].[68] F. Bahr, D. Banerjee, F. Bernardoni, M. Koren, H. Simma and R. Sommer, Extraction of bareform factors for B s → K ‘ν decays in nonperturbative HQET , Int. J. Mod. Phys. A (2019)1950166, [ ].[69] H. Akaike, Information theory and an extension of the maximum likelihood principle , in
Springer Series in Statistics , p. 199. Springer New York, 1998.[ doi:10.1007/978-1-4612-1694-0_15 ].[70] J. Heitger and A. Jüttner,
Lattice cutoff effects for F D s with improved Wilson fermions - afinal lesson from the quenched case , JHEP (2009) 101, [ ]. [Erratum: JHEP (2011) 036].[71] W. I. Jay and E. T. Neil, Bayesian model averaging for analysis of lattice field theory results , .[72] G. Schwarz, Estimating the dimension of a model , Ann. Statist. (03, 1978) 461.[73] C. M. Hurvich and C.-L. Tsai, Regression and time series model selection in small samples , Biometrika (06, 1989) 297.[74] M. Dalla Brida, R. Höllwieser, F. Knechtli, T. Korzec, A. Ramos and R. Sommer, Non-perturbative renormalization by decoupling , Phys. Lett. B (2020) 135571,[ ].[75] G. S. Bali, S. Collins, F. Hutzler, M. Göckeler, A. Schäfer, E. E. Scholz et al.,
Towards thecontinuum limit with improved Wilson fermions employing open boundary conditions , PoS
LATTICE2016 (2017) 106, [ ].
76] G. S. Bali et al., in preparation.[77] M. Bruno, M. Dalla Brida, P. Fritzsch, T. Korzec, A. Ramos, S. Schaefer et al.,
QCDCoupling from a Nonperturbative Determination of the Three-Flavor Λ Parameter , Phys. Rev.Lett. (2017) 102001, [ ].[78] D. Hatton, C. Davies, B. Galloway, J. Koponen, G. Lepage and A. Lytle,
Charmoniumproperties from lattice QCD + QED : Hyperfine splitting,
J/ψ leptonic width, charm quarkmass, and a cµ , Phys. Rev. D (2020) 054511, [ ].[79] F. Bernardoni et al.,
The b-quark mass from non-perturbative N f = 2 Heavy Quark EffectiveTheory at O (1 /m h ), Phys. Lett. B (2014) 171, [ ].[80] K. Chetyrkin, J. H. Kühn and M. Steinhauser,
RunDec: A Mathematica package for runningand decoupling of the strong coupling and quark masses , Comput. Phys. Commun. (2000)43, [ hep-ph/0004189 ].[81] B. Schmidt and M. Steinhauser,
CRunDec: a C++ package for running and decoupling of thestrong coupling and quark masses , Comput. Phys. Commun. (2012) 1845, [ ].[82] F. Herren and M. Steinhauser,
Version 3 of RunDec and CRunDec , Comput. Phys. Commun. (2018) 333, [ ].[83] (Particle Data Group) , P. Zyla et al.,
Review of Particle Physics , PTEP (2020)083C01.(2020)083C01.