Continuum limit of the leading order HQET form factor in $B_s \to K\ell\nu$ decays

We discuss the computation of form factors for semi-leptonic decays of $\rm B$-, $\rm B_s$- mesons in lattice QCD. Considering in particular the example of the static $\rm B_s$ form factors we demonstrate that after non-perturbative renormalization the continuum limit can be taken with confidence. The resulting precision is of interest for extractions of $V_{\rm ub}$. The size of the corrections of order $1/m_{\rm b}$ is just estimated at present but it is expected that their inclusion does not pose significant difficulties.


Introduction
Weak decays of B-mesons are a very important piece in the puzzle of understanding about how well the Standard Model of particle physics describes Nature. One relevant question concerns the determination of the Cabibbo-Kobayashi-Maskawa matrix element V ub from different decays. This fundamental parameter of the Standard Model is not known very precisely yet. Testing for consistent values provides a check of the Standard Model. In fact results extracted from inclusive decays agree with those from different exclusive decays, like B → π ν or B → τ ν [1][2][3], only after stretching the presently estimated uncertainties by around a factor three. We avoid calling this a three-sigma tension since the uncertainties are largely systematic, coming from the theoretical computation of form factors in lattice QCD on one side and the perturbative treatment of inclusive decays on the other side. But also experimental uncertainties contribute.
In this letter we consider the determinations of semi-leptonic form factors for B s -mesons from lattice QCD. A review with some discussion of the challenges involved is found in [4]. It appears that the most relevant challenge is the presence of a (large) mass scale m b ∼ 5 GeV. Together with inverse lattice spacings below 4 GeV this distorts the continuum physics considerably in a straight application of lattice QCD. We do not want to review here this issue in detail, but just mention that this leads one to consider effective field theories for the b-quark or extrapolations in its mass, again guided by effective field theory considerations. The most advanced computations [5][6][7][8], use either a relativistic heavy quark action or employ non-relativistic QCD on the lattice. There the challenge is twofold. First, a fully non-perturbative renormalization program for the heavy-light currents does not (yet) exist. It is replaced by "mostly non-perturbative" renormalization [9,10], where the factor Z hl / √ Z hh Z ll is taken from 1-loop perturbation theory and this approximation is expected to be a good one [9,10]; alternatively straight 1-loop perturbation theory is used. Second, discretisation errors are estimated only by power-counting arguments because continuum limit extrapolations may involve a complicated functional dependence on the lattice spacing. As a consequence we are not aware of the computation of a non-perturbatively renormalized heavy-light form factor extrapolated to the continuum.
In order to place the present work into context, let us briefly list the steps which are necessary to come to a trustworthy result of interest to phenomenology: a) obtain the ground state matrix elements K|V µ (0)|B s that mediate the transition, b) renormalize the currents (and thus matrix elements) and, if an effective theory is used, relate them to QCD ("matching"), c) take their continuum limit, d) extrapolate to the quark masses realized in Nature, e) map out the q 2 dependence.
Here we demonstrate our solutions to a)-c). We are very brief about our specific choices in a), even though the extraction of the ground-state-to-ground-state matrix element which gives the leading form factor, h ⊥ , is delicate since excited state contributions have significant amplitudes. Details on this more technical but important issue are relegated to a companion paper [11]. Steps d)-e) will follow in due course.
We concentrate on the non-perturbative renormalization and the continuum limit albeit only in the leading order of Heavy Quark Effective Theory (HQET). We consider a single value of the momentum transfer and a single value of the light dynamical quark masses with only two degenerate dynamical flavors. These restrictions mean that our computation does not immediately advance phenomenology, but since a continuum limit was not taken before we can see for the first time how it works, and the result provides a cross check on the uncertainty estimates of previous computations. We will find that with our discretisation and with non-perturbative renormalization, the continuum limit (for Kaon momentum around 0.5 GeV) is smooth. Given that the inclusion of 1/m b effects in the systematic treatment of HQET [12] was not a severe problem (apart from a lot of work) in other quantities [13][14][15] we are very encouraged to complete the started programme towards phenomenologically relevant results.

Form Factors
We consider the decay B s → K ν. Working at the leading order in the weak interactions, the transition amplitude factorises into a straightforward leptonic amplitude and the QCD matrix element with the equivalent form-factor decompositions The last line, with defines the form factors, h and h ⊥ . The usual squared momentum transfer q 2 = (q 0 ) 2 −q 2 , with q µ ≡ p µ Bs − p µ K , is here replaced by which at fixed Kaon four-momentum p K is independent of the mass of the b-quark. This property, together with the factor (2m Bs ) 1/2 in (2.1), is convenient to discuss the behavior of the amplitude at large mass of the b-quark. It removes the mass-dependence of the (standard) relativistic normalization B s (p )|B s (p) = 2E(p)(2π) 3 δ(p − p ) of the state of the heavy meson. Since the current V µ (x) ≡ψ u (x)γ µ ψ b (x) translates into heavy quark effective (mass-independent) fields with only a logarithmically mass-dependent conversion function, the form factors h and h ⊥ scale only logarithmically with the mass in the limit of large b-quark mass.
Choosing for the remainder of this letter the rest-frame of the B s -meson with v µ = (1, 0, 0, 0) as a reference frame, the invariant kinematic variable is just p K · v = E K , the energy of the final-state pseudo-scalar. Upon neglecting terms proportional to m 2 /m 2 Bs and m 2 /q 2 (m being the mass of the final-state lepton), the differential decay rate is then given by A comparison of (2.3) with the experimentally measured rate allows for a determination of |V ub | once the form factors are known. They need to be determined at a single value (or ideally in a range) of E K where overlap with experimental data exists. In our frame (p Bs = 0), the form factors are obtained from the (QCD) matrix elements With the above normalization, eq. (2.1), they have an HQET expansion without factors that involve a power of the quark mass. Rather the r.h.s. depend logarithmically on the mass of the heavy quark, due to the matching of HQET to QCD. In our notation, are the renormalization group invariant (RGI) operators in HQET, and the conversion functions C x connect (the matrix elements of) V stat,RGI 0 and V stat,RGI k to the ones in QCD, see [16,17].
The functions C x are known with 2-loop matching (for short "2-loop"), i.e. up to α(m b ) 3 corrections in continuum perturbation theory [18][19][20][21][22][23][24][25]. We use them here with the RGI b-quark mass M b and the Λ-parameter determined in the theory with two dynamical flavors [13,26], i.e. M b /Λ MS = 21.2(1.2), where the uncertainty of Λ dominates. The conversion functions then evaluate to 1 where the second quoted uncertainty is estimated as the difference between 2-loop and 1-loop. It is not entirely clear whether this is a conservative estimate of the perturbative error (see sect. 2.3 of [12]). Let us clarify the difference to standard 1-loop renormalization of heavy-light form factors. We renormalize the HQET currents in (2.8) non-perturbatively, thus the continuum limit of their matrix elements is not affected by any perturbative uncertainty. Then, in the continuum, the factor C x is known only perturbatively, but to one more power of α s than what is available for the total Z x = C x × Z stat,RGI x in other approaches. Thus, even if we quote an uncertainty of up to five percent for the renormalization, this is an O(α 3 ) uncertainty where usually it is O(α 2 ). In the future the ALPHA collaboration will non-perturbatively match HQET and QCD [28] also for the vector currents [29]. One then obtains directly Z x = C x × Z stat,RGI x with full non-perturbative precision. We now proceed to the numerical evaluation of the m b -independent RGI matrix elements h stat,RGI x , which are not affected by perturbative errors or ambiguities. 1 They are conveniently summarized in [27].

Framework and Renormalization
For our first numerical investigation of the problem, we choose N f = 2 flavors of quarks. The prime reason for this choice is that in a related project, the non-perturbative matching of HQET to QCD for the currents V 0 , V k is being carried out at the order 1/m b [29][30][31][32]. Once this is complete, we will be able to include the 1/m b corrections with little additional effort. For now, we remain at the lowest order of HQET, namely the static order.
The b-quark is then replaced by a static quark [33] labelled "h". Two different discretisations, HYP1 and HYP2, are chosen [34]. These have moderate discretisation errors and a much improved signal-to-noise ratio compared to the classic Eichten-Hill static quark action [33]. The bare currents are form-identical to the ones in QCD, apart from the O(a) improvement terms ( ← − ∇ S denotes the symmetric covariant derivative acting on the field to the left). The coefficients 2 , V 0 = 0.0380 (6) .
The multiplicative renormalization of the currents eq. (2.8) can be written as where we have made explicit that Z stat,RGI renormalizes exactly as A stat 0 due to the spin symmetry of HQET, while for V 0 an extra factor Z stat V/A (g 0 ) originates from the broken chiral symmetry of Wilson fermions. We use the non-perturbative results relating A stat 0 , renormalized at scale µ 0 = 1/L max in the Schrödinger functional scheme, to the RGI operator. The factor R is known non-perturbatively by running to very high µ and continuum extrapolation. In Table 1, the remaining piece Z stat A (g 0 , aµ 0 ) is reproduced from [36]. For the finite renormalization Z stat V/A we use a range 2 Spin symmetry leads to the identity cV k = cA 0 . Ref. [35] uses the notation c stat ( 5) [38], an update of [26]. All ensembles have T = 2L and m π L ≥ 4 where m π is the pion mass. N cfg denotes the number of configurations on which we performed measurements. The hopping parameter of the strange quark is denoted by κ s . The angle θ appears in the flavor-twisted boundary conditions, as explained in the text.
This range is generous, because Z stat V/A has been seen to be very close to one in the quenched approximation [37], and at 1-loop order, there is no N f -dependence. Despite these arguments our range in (3.7) is no more than an educated guess. This is adequate since Z stat V/A affects only one of the 1/m b suppressed terms, which we just use as an illustration of the associated uncertainty. Again, we note that such issues will be eliminated when the non-perturbative matching is carried out.

Simulation Parameters
We base our investigation on a subset of the lattice gauge field configurations with two degenerate flavors of improved Wilson fermions and Wilson gauge action generated by the Coordinated Lattice Simulations (CLS) effort [26]. The observables are computed on three ensembles, namely A5, F6 and N6, chosen to have roughly the same pion mass but three different lattice spacings, see Table 2. The quoted lattice spacings were determined from the Kaon decay constant f K in reference [26] and updated in reference [38].
Our choice of the gauge field ensembles fixes the degenerate up and down quark masses, m l . For the spectator strange (valence) quark mass we are free, however, to choose any smooth function m s (m l ) which passes through the physical point. As in [26] we define this function by fixing the Kaon mass in units of the Kaon decay constant to its physical value at our (and any) value of m l . We expect that this will lead to a flat extrapolation to the physical value of m l , the "physical point".
We choose |p K | = 0.535 GeV which corresponds to the minimum available momentum with periodic boundary conditions for all fields on the N6 lattice. On the other lattices we keep p K = (1, 0, 0) (2π + θ)/L fixed by introducing flavor-twisted boundary conditions  Figure 1: Effective energy of C K (left) and C Bs (right) on ensemble N6. Note that both panels have equal ranges on the corresponding axes. One can identify reasonable plateaus with small errors even though the Kaon carries a non-vanishing momentum and we have a static B s -meson. The value for the ground state energy as obtained from a two-exponential fit is shown as a red band. Uncertainties shown here are only those of E stat and E K in lattice units, not the ones of m bare and the lattice spacing. The data points are shown for the case of maximum smearing, while the fit involves all the smearing levels.
[39] ψ s (x + L1) = e iθ ψ s (x), for the strange quark. The B s -meson is kept at rest by ψ h (x + L1) = e iθ ψ h (x), and we remain with periodic boundary conditions for all other fields. The numerical values for θ are listed in Table 2. Our choice of p K yields a central value of q 2 = 21.22 GeV 2 at all lattice spacings and an error coming from the lattice spacing of 0.03 − 0.05 GeV 2 . Note that the flavor-twist is introduced only for quenched quarks.

Correlation Functions and Matrix Elements
We work with all-to-all light quark propagators [40,41] implemented by a random U(1) source placed on each time slice ("full time dilution" [42]). While this is numerically costly, it significantly reduces the large-time variance. Together with the deflated solver [43][44][45] that we use, it is thus very cost effective. Most notably, this feature of our computation provides access to all time separations of two-point and three-point functions on our lattices. For details we refer to [11].
The two-point functions are defined as with the pseudoscalar density P (i) . Indices i and j denote the use of different Gaussian wave functions for the light quarks with a set of smearing parameters as in [13]. Up to terms which are exponentially suppressed in the time extent T of the torus, we can parameterise C K , C Bs as where we have denoted amplitudes by κ (n) = L 3/2 (2E K ) −1/2 0|P us (0)|K, n and β Bs . The restrictions t x a are required because we use an improved action where the positivity of the transfer matrix is not guaranteed. In a fit to these correlation functions, we use t-ranges denoted by t K2 min ≤ t K ≤ t K2 max and t B2 min ≤ t Bs ≤ t B2 max , respectively, and N is the number of excited B s -meson states which we include. Note that we have restricted ourselves here to only the ground state of the Kaon. There is a single, fixed, smearing level for the Kaon. In the Kaon two-point function, ground-state dominance sets in at around 1.2 fm, while for the B s -meson this happens a bit earlier for our optimal (widest) Gaussian wave function. An illustration is found in Figure 1, which also shows that a reasonably good precision is reached; plateaus are also clearly visible at the larger lattice spacings.
The three-point function has a representation C Bs→K We perform a simultaneous fit to C K , C Bs ij , and C Bs→K µ,i . For the latter we consider (t K , t Bs )values restricted to the rectangle with t K3 min ≤ t K ≤ t K3 max,µ and t B3 min ≤ t Bs ≤ t B3 max,µ . The desired form factors are given by the ground-state matrix elements Their extraction from the data is rather delicate because statistical errors grow with time separations. Due to our all-to-all computation and the use of HYP1, HYP2 discretisations, we still have a precision of better than two percent for t K 2 fm and t Bs 1.5 fm in the two-point functions; for the three-point function it drops below the two-percent level at around t K + t Bs ≈ 2 fm (in the fitting region). However, we find that N = 2 excited states are necessary in our fits to obtain a good description of the data and, in particular, a safe extraction of the most important form factor ϕ (0,0) k . To determine the reliability of the fits, we vary the boundaries t B2 min , t K3 min , and t B3 min of the fit ranges and verify that the fit results remain unchanged within errors. As an example, Fig. 2 shows the dependence of the fit results for ϕ (0,0) 0 and ϕ (0,0) k on t K3 min and t B3 min (keeping t B3 min − t B2 min = 5a fixed). The other boundaries are chosen to suppress the effects of excited Kaon states (t K2 min ), of noise (t B2 max ), and of the finite time extent T . The latter two criteria considerably constrain our choice of t K3 max,µ and t B3 max,µ . Table 3 lists the fit ranges which we used for the HYP2 data. The bare ground-state matrix elements are shown in Table 4. Further details will be described in [11].

Continuum Limits
The bare form factors, renormalized as explained in Sect. 3.1, yield the RGI form factors listed in Table 5. Their errors take account of all statistical correlations and autocorrelations as described in [13] based on [46,47]. For these non-perturbatively renormalized form factors (at fixed squared momentum transfer q 2 , or Kaon energy p K · v) the continuum limit can now be taken. Fig. 3 shows the dependence of the results on the lattice spacing and the continuum extrapolation for the two discretisations (HYP1 and HYP2). As our best result we estimate the continuum limit by a linear extrapolation in a 2 of the data with c x = c (1) x g 2 0 , as illustrated in the figure.
For h a simple constant extrapolation (a weighted average) yields compatible results, but of course with much smaller error bars. Since there is no reason, why a 2 effects should be entirely absent, we use the numbers with the larger error bars.
It seems not critical that we know the O(a) improvement coefficients of the currents only in 1-loop perturbation theory. These coefficients are not very relevant at the level of precision of our data. 3 Finally we combine the continuum limits of HYP1 and HYP2 in a weighted average,

The Form Factors f + and f 0 , and a Comparison to Other Results
We now have different options to estimate the form factor f + . Working only at static order, we can use any quantityf +,i with and we consider (dropping the argument of where the second errors are the ones from the perturbative uncertainty in (2.9) and (2.10).
In order to give a single number, we estimate the O(1/m b ) terms in this way and quote f + (21.22 GeV 2 ) =f +,2 ± 0.15f +,2 = 1.63(8)(6) ± 0.24 . Besides f + it is common in the literature to report results for the scalar form factor f 0 , which is another linear combination of h ⊥ and h . To estimate its value, we usẽ Our results, eq. (3.19) and eq. (3.21), compare rather well with other values of the form factors in the literature. The result of Flynn et al. [6], extracted at our value of q 2 , is f + ≈ 1.65 (10) and f 0 ≈ 0.62(5), while Bouchard et al. [7] have f + ≈ 1.80 (20) and f 0 ≈ 0.66 (5). As should be clear, our estimates have a systematic error of a completely different nature. While we focused our effort on taking the continuum limit of non-perturbatively renormalized matrix elements at a fixed Kaon momentum, we have so far neglected the dependence on the light-quark mass which -according to [6,7] -is below our errors.

Conclusion and Outlook
For the first time we have been able to perform a study of the continuum limit of fully non-perturbatively renormalized form factors. They are computed at a fixed squared momentum transfer q 2 = 21.22 GeV 2 , and we have concentrated on the leading-order form factors in 1/m b . In RGI form these are unambiguous. Our main result is contained in Fig. 3. It shows that the continuum limit at a Kaon momentum of 1/2 GeV is smooth. This behavior of the discretisation errors for matrix elements with a momentum of this size is not obvious a priori and linked to our choice of actions; a generalization is at most possible at a rather qualitative level. With this encouraging result, we may consider also somewhat larger momenta in the future.
Our numbers in eq. (3.19) and eq. (3.21) provide a positive cross-check of [6,7]. Already now, they thus strengthen our confidence in the form factors extracted on the lattice and summarized in [4], but they will be of a more direct phenomenological interest when the 1/m b terms are included and the errors shrink accordingly. At that point we also have to carefully consider the extrapolation to the physical light quark masses and finally more than one value of the Kaon momentum will be of interest. As a bottom line, the study presented here suggests that all this is possible with a precision which is of interest for an extraction of V ub from experimental decay rates.