The b-quark mass from non-perturbative $N_f=2$ Heavy Quark Effective Theory at $O(1/m_h)$

We report our final estimate of the b-quark mass from $N_f=2$ lattice QCD simulations using Heavy Quark Effective Theory non-perturbatively matched to QCD at $O(1/m_h)$. Treating systematic and statistical errors in a conservative manner, we obtain $\overline{m}_{\rm b}^{\overline{\rm MS}}(2 {\rm GeV})=4.88(15)$ GeV after an extrapolation to the physical point.


Introduction
The masses of the quarks are among the fundamental parameters of the Standard Model (SM), and as such hold considerable interest. Heavy quark masses, in particular, enter as parameters in various perturbative predictions of interesting decay rates, e.g.
B → X s γ or inclusive B → u or B → c rates. Such decays yield useful constraints for the CKM matrix and, in principle, options to obtain hints for physics beyond the Standard Model. It is therefore desirable to minimize the uncertainty in m b entering these predictions.
The b-quark mass also enters the prediction for the cross section of the H → bb decay, which is the mode with the largest branching ratio for an SM-like Higgs with a mass of 126 GeV. In the future, tests of this coupling will help providing further characterizations of the new boson.
lattice approach in general offers the unique opportunity to study the N f -dependence of the b-quark mass. We will further discuss that in the conclusions.
In this Letter, we present our results for the mass of the b-quark from simulations of non-perturbative N f = 2 HQET. In Section 2, we briefly review the methods employed before presenting the results in Section 3. Section 4 contains our conclusions.

Methodological background
HQET on the lattice constitutes a theoretically sound approach to heavy quark physics by expanding QCD correlation functions into power series in 1/m h around the static limit m h → ∞, which is non-perturbatively renormalizable, so that the continuum limit can always be taken.
Following the strategy described in [35,39] and previously applied to calculate m b in the quenched approximation [36], we write the HQET action at O(1/m h ) as (2.4) where the heavy quark spinor field ψ h obeys 1+γ 0 2 ψ h = ψ h , m bare is the bare heavy quark mass absorbing the power-divergences of the self-energy in the static approximation, and the parameters ω kin and ω spin are formally of order 1/m h and have been previously determined in [38].
for the expectation value of some multilocal fields O , where O stat is the expectation value of O determined in the static theory. A significant improvement in the signal-to-noise ratio of heavylight correlation functions can be achieved by defining the co- )/a in terms of a suitably smeared link instead of the bare link U (x, 0). Each smearing prescription constitutes a separate lattice action; here we have employed both the HYP1 and HYP2 actions [40][41][42].
To reliably extract hadronic quantities, we have to pay particular attention to unwanted contributions from excited states. The variational method, which has become a standard tool for analyzing hadronic spectra in lattice QCD, starts from correlator matrices C stat The main ingredient is to solve the generalized eigenvalue problem (GEVP) in the static limit Indeed, by exploiting the orthogonality property of the eigenvectors v stat one can show that the O(1/m h ) corrections to the energy levels depend only on the static generalized eigenvalues λ stat n (t, t 0 ), the eigenvectors v stat n (t, t 0 ), and the O(1/m h ) correlators C kin/spin (t) [43], in analogy with perturbation theory in quantum mechanics. At large times t and t 0 satisfying t 0 t/2, the asymptotic behaviour is then known to be [43] E eff,stat (2.9) with E m,n = E m − E n . The time intervals over which we fit the energy plateaux are chosen so as to minimize the systematic error from the excited states while keeping the statistical error under control.
Finally, the mass of the B-meson to O(1/m h ) is given by

Ensembles used
Our measurements are carried out on a subset of the CLS (Coordinated Lattice Simulations) ensembles, which have been generated using either the DD-HMC [44][45][46][47] or the MP-HMC [48] algorithm, using the Wilson plaquette action [49] and N f = 2 flavours of non-perturbatively O(a) improved Wilson quarks [50,51]. An overview of the simulation parameters of the ensembles used is given in Table 1. In order to suppress finite-size effects, we consider only ensembles satisfying m π L > 4.0. The light valence quarks are equal to the sea quarks, and the (quenched) b-quark is treated by HQET.
In order to control the statistical error in a reliable fashion, we make use of the method of [53] as improved by [52] to estimate the effect of long-term autocorrelations due to the coupling of our observables to the slow modes of the Markov chain, decaying as ∼ exp(−τ /τ exp ) in Monte Carlo simulation time τ . The propagation of these effects through to the continuum-extrapolated result at the physical pion mass is carried out iterating the formulae of [52].

Lattice spacings
The lattice spacings a, pion masses m π and pion decay constants f π on the CLS ensembles used here are taken from an update [54] of the analysis in [55] with increased statistics and including additional ensembles. They read a = 0.04831 (38) Table 1 Details of the CLS ensembles used: bare coupling β = 6/g 2 0 , lattice spacing a, spatial extent L in lattice units (T = 2L), pion mass m π , m π L, number of configurations employed, and number of configurations employed normalized in units of the exponential autocorrelation time τ exp as estimated in [52]. Additionally, we specify the CLS label id and the Gaussian smearing parameters R k used to build different interpolating fields as described in the text.  that is needed to convert the b-quark mass into physical units later on. The length scale L 1 originates from the non-perturbative finite-volume matching step used to determine the HQET parameters [38].

Basis of B-meson interpolating fields
Our basis of N = 3 operators is given by where ψ h (x) is the static quark field, and different levels of Gaussian smearing [56] with a triply (spatially) APE smeared [57,58] covariant Laplacian are applied to the relativistic quark field ψ (k) Our smearing parameters κ G = 0.1 and R k , collected in Table 1 3 fm, and t max is fixed to ∼ 0.9 fm. This will assure that our selection criterion σ sys σ /3 is satisfied [59], is shown in Fig. 1.

Determination of the b-quark mass
The mass of the B-meson to static order is given by  Table 2.
Once m B,δ have been computed for a set of z spanning a range of heavy quark masses containing the b-quark mass, we perform a combined chiral and continuum extrapolation to obtain m B (z, m exp π ) ≡ m B,δ (z, m exp π , 0), using m exp π = 134.98 MeV [1].
Considering that the O(a) improvement was performed nonperturbatively but neglecting O(a/m b ) effects, 1 the NLO formula from HMChPT reads [60] m sub the leading non-analytic term of HMChPT has been subtracted. The B * Bπ couplingĝ = 0.489 (32) has been determined recently [61] and the variable y is identical toỹ 1 introduced in [55]. We use the convention where the pion decay constant is f exp π = 130.4 MeV.
The extrapolation (3.7) is shown in Fig. 2 (left) for three values of z in the vicinity of z b = M b L 1 we are aiming at. Its result, B(z) = m B,δ (z, m exp π , 0), is given in Table 2. As shown in Fig. 2 1 Accounting for an a/m b has little effect. Adding a term F δ · (a/m b ) to Eq.
(3.12) 2 We follow the notation of Gasser and Leutwyler [62] for the definition of the with a conversion function ρ(r) that can be evaluated accurately using the known 4-loop anomalous dimensions of quark masses and coupling [63,64]. It is described in more detail in Appendix A.
The ratio r b = M b /Λ MS is computed from our value of z b and the ALPHA Collaboration results for non-perturbative quark mass renormalization [65]. We find r b = 21.1(13), ρ(r b ) = 0.640 (6)  We emphasize that this is the mass in the theory with two dynamical quark-flavours, the b-quark is quenched, a completely well defined approximation for a heavy quark. In particular also the function ρ(r) refers to N f = 2.
To have an idea of the magnitude of O(1/m b ) corrections to the b-quark mass one must repeat the above computation in the static limit. The reason is that the 1/m contribution ω kin E kin + ω spin E spin is divergent in the continuum limit; only the combination with m bare in Eq. (2.10) is finite. The HQET parameter m stat bare was determined by matching the static theory with QCD as described in [38]. By repeating the same steps as for the NLO case we ob-  The result of the combined chiral and continuum extrapolation of m B in the static limit, as well as the quadratic interpolation Table 3 Partial contributions (σ i /σ ) 2 to the accumulated error σ of z b . Only error sources contributing with a relative squared uncertainty (σ i /σ ) 2 > 0.5% are listed. The ensemble A3 did not appear in Table 1 since it enters through the scale setting procedure [54,55] only.
1.2 0.9 2.6 5.9 5.6 20. 6 61.6 in z to obtain z stat b are very similar to those obtained at next-toleading order. The small differences observed between the results in (3.10)-(3.11) and (3.14) show that for this observable the HQET expansion is very precise, making us confident that O(1/m 2 b ) corrections are negligible with present accuracy. Indeed, the smallness of the 1/m b terms is known with much higher accuracy than (3.10) suggests, e.g., z We conclude this section by analyzing the error budget for z b . As can be seen in Table 3 approximately 62% of the contribution to the square of the error comes from the HQET parameters.
Another ∼ 21% comes from the relativistic Z A that affects the computation of z b through the scale setting, while only the residual ∼ 17% comes from the computation of the HQET matrix elements.
In this respect the largest contribution comes from the ensembles at β = 5.5, that are more affected by long-term autocorrelations (critical slowing down).

Discussion and conclusions
Using non-perturbatively matched and renormalized HQET in N f = 2 lattice QCD, we have determined the mass of the b-quark with essentially controlled systematic errors: in particular, the renormalization is carried out without recourse to perturbation theory and the continuum limit is taken. An irreducible systematic error which remains is a m b /m b ∼ (Λ/m b ) 3 relative error due to the truncation of the HQET expansion at order Λ 2 /m b . However, with a typical scale of Λ = 500 MeV one obtains a permille-sized truncation error, which is completely negligible with today's accuracy. The estimate is supported by the fact that we do not see any difference between our static result and the one including the Λ 2 /m b terms. Furthermore, according to previous experience an effective scale of around Λ = 500 MeV seems to govern the expansion [38,59,66].
Our results, are in agreement with the N f = 2 results of [28] who cite a similar error, but use a completely different approach. We compare to the quenched approximation and to the PDG values in Table 4. There is little dependence of m MS b (μ) on the number of flavours for N f = 0, 2, 5 and for typical values of μ between m MS b itself and 2 GeV. In particular at the lower scale of 2 GeV, where the apparent convergence of perturbation theory is still quite good, a flavour number dependence of the mass of the b-quark is not detectable at all. In hindsight, this is rather plausible as we match our effective theories (albeit with only N f = 0, 2 dynamical flavours) to the real world data at low energies. Indeed, precisely speaking the above statements refer to the theories renormalized by fixing the B-meson mass to its physical value and setting the overall energy scale through the kaon decay constant [55] or roughly equivalent the pion decay constant [54]. 3 In this way the low energy hadron sector of the theories is matched to experiment, and it is natural to expect that the quark masses agree at a relatively low scale. On the other hand we do not want to push the perturbation theory needed for giving m b in the MS scheme to scales below 2 GeV. We remark that also the strange quark mass at 2 GeV is known to be only weakly dependent on N f [67,68].
In contrast, the RGI mass M b differs significantly between N f = 5 and N f = 2. Given the observed weak flavour number dependence at scales of 2-5 GeV, the differences in M b can be traced back to the N f dependence of both the RG functions and the Λ parameters. These two effects happen to reinforce each other between N f = 5 and N f = 2 while in the comparison N f = 2 and N f = 0 they partially compensate.
All of this suggests to use the b-quark mass at scales around μ = 2 GeV when one attempts to make predictions from theories with a smaller number of flavours for the physical 5-flavour theory. With a less detailed look, the overall picture of the MS masses in Table 4 suggests that -at the present level of errors -the b-quark mass is correctly determined from the different approaches. Our method is very different from those which enter the PDG average. It avoids perturbative errors in all stages of the computation except for the connection of the RGI mass to the running mass in the MS scheme, where truncation errors seem to be very small. Due to these properties, it remains of interest to apply our method with at least three light dynamical quarks and test the consistency of the table once more. As remarked earlier, the error budget of our present computation is such that in a future computation a significantly more precise number can be expected. Table 4 Masses of the b-quark in GeV in theories with different quark flavour numbers N f and for different schemes/scales as well as Λ MS and the RGI mass M. The PDG value of the b-quark mass is dominated by [8,24].  termined by the running of the quark mass and the coupling in the Schrödinger Functional (SF) scheme. They are known nonperturbatively in terms of the step scaling functions of [65]. For the error analysis we take the errors in h, k including their correlation into account, remembering that h also contributes through . The uncertainty arising from the perturbative running in the MS scheme is negligible. For example, adding the recently computed 5-loop term in the mass anomalous dimension [71] does not change numbers at the one permille level. The error analysis for m MS (μ) with some fixed μ is carried through analogously.