'Bs -->Ds l nu' near zero recoil in and beyond the Standard Model

We compute the normalization of the form factor entering the Bs -->Ds l nu decay amplitude by using numerical simulations of QCD on the lattice. From our study with Nf=2 dynamical light quarks, and by employing the maximally twisted Wilson quark action, we obtain in the continuum limit G(1) = 1.052(46). We also compute the scalar and tensor form factors in the region near zero recoil and find f0(t0)/f+(t0)=0.77(2), fT(t0,mb)/f+(t0)=1.08(7), for t0=11.5 GeV^2. These latter results are useful for searching the effects of physics beyond the Standard Model in Bs -->Ds l nu decays. Our results for the similar form factors relevant to the non-strange case indicate that the method employed here can be used to achieve the precision determination of the B -->D l nu decay amplitude as well.

We compute the normalization of the form factor entering the Bs → Ds ν decay amplitude by using numerical simulations of QCD on the lattice. From our study with N f = 2 dynamical light quarks, and by employing the maximally twisted Wilson quark action, we obtain in the continuum limit G(1) = 1.052(46). We also compute the scalar and tensor form factors in the region near zero recoil and find f0(q 2 0 )/f+(q 2 0 ) = 0.77 (2), fT (q 2 0 , m b )/f+(q 2 0 ) = 1.08 (7), for q 2 0 = 11.5 GeV 2 . These latter results are useful for searching the effects of physics beyond the Standard Model in Bs → Ds ν decays. Our results for the similar form factors relevant to the non-strange case indicate that the method employed here can be used to achieve the precision determination of the B → D ν decay amplitude as well.

I. INTRODUCTION
Inclusive and exclusive semileptonic b → c ν decays, with ∈ {e, µ}, have been subjects of intensive research over the past two decades. Within the Standard Model (SM) the main target of that research was, and still is, the accurate determination of the Cabibbo-Kobayashi-Maskawa matrix element |V cb |, which is extracted from the comparison of theoretical expressions with experimental measurements of the partial or total decay widths. It turns out, however, that the value for |V cb | obtained from the exclusive decays agrees only at the 2.1 σ level with the one extracted from inclusive decays. More specifically, the two independent analyses of inclusive decays [1,2] (updated in refs. [3]) give completely consistent results which, after taking an average, reads 1 |V cb | incl. = 41.90(70) · 10 −3 . (1) The analyses of exclusive decays, instead, are performed by fitting the experimental data to the shapes of the form factors parameterized according to the expressions proposed and derived in ref. [5], so that the final results are then reported in the following form, as obtained from B(B → D * µν) and B(B → Dµν), respectively. F(1) and G(1) are the relevant hadronic form factors at the zero recoil point. Thanks to the heavy quark symmetry, and up to perturbative QCD corrections, both these form factors are equal to one in the limit of m c,b → ∞ [6].
To compute the deviation of these form factors 1 The most recent update of the analysis of inclusive decays using the kinetic scheme can be found in ref. [4].
from unity, it is necessary to include all the nonperturbative order 1/m n c,b QCD corrections. The only model independent method allowing to compute F(1) and G(1) from the first theory principles is lattice QCD. Using the most recent estimates of the above form factors, F(1) = 0.902 (17) [7] and G(1) = 1.074 (24) [8], one arrives to |V cb | excl. = 38.56(89) · 10 −3 , a number that obviously differs from |V cb | incl. given in eq. (1). In principle the exclusive decay modes are better suited for the precision determination of |V cb | because less theoretical assumptions are needed to compute the corresponding decay rates. The main obstacle is, however, the necessity for a reliable, high precision, lattice QCD estimate of F(1) and G (1). Furthermore, for the required percent precision of |V cb | it is important to have a good control over the structure dependent soft photon B → D ( * ) νγ soft which could otherwise be misidentified as pure semileptonic decays. This problem is less acute for the down-type spectator, i.e.B 0 → D ( * )+ ν, than for the up-type spectator quark, B − → D ( * )0 ν [9]. For that reason it is desirable to consider the charged and neutral Bmeson decay modes separately. In this paper we will mainly discuss B s → D s ν decay for which the soft photon pollution is smaller. Moreover, this mode is much more affordable numerically because the valence s-quark is easily accessible in numerical simulations of QCD on the lattice which is not the case with the physical u/d-quark. In this paper we will also comment on the non-strange case when appropriate. Finally, we prefer to focus on B (s) → D (s) ν, rather than B (s) → D * (s) ν, because the hadronic matrix element involves much less form factors and the decay rate is therefore likely to be less prone to systematic uncertainties.
Despite the importance of G(1), only a few lattice QCD studies have been performed so far. The methodology described in ref. [10] has been implemented in unquenched simulations with N f = 2 + 1 dynamical staggered light quarks in refs. [8,11] where the propagating heavy quarks on the lattice have been interpreted by means of an effective theory approach. In ref. [12], an alternative method to compute the B → D semileptonic form factors has been proposed and implemented in quenched approximation. That latter method, based on use of the step scaling function, allows to compute the same form factors without recourse to heavy quark effective theory. The price to pay, however, is that the method of ref. [12] is computationally very costly and to this date it has not been extended to unquenched QCD. In this paper we use a mod-ification of the proposal of ref. [12], presented in ref. [13] and also implemented in the computation of the decay constant f B and the b-quark mass, cf. ref. [14]. The remainder of this paper is organized as follows: In sec. II we define the form factors and express them in a way that is suitable for the strategy used for their computation which is described in sec. III. Details of our lattice computations and the results for G(1) are given in sec. IV, while the results concerning the scalar and tensor form factors in the region close to zero recoil are discussed and presented in sec. V. We finally conclude in sec. VI.

II. DEFINITIONS
The hadronic matrix element describing the B s → D s ν decay in the SM, D s |bγ µ (1 − γ 5 )c|B s ≡ D s |bγ µ c|B s , is parameterized in terms of the hadronic form factors f +,0 (q 2 ) as where V µ =bγ µ c, q = p − k, and q 2 ∈ (0, q 2 max ], with q 2 max = (m Bs − m Ds ) 2 . The extraction of the form factors becomes particularly simple if we use the projectors so that In our computations we will consider the situations with | p| = 0 and q µ = (m Bs − E Ds , − k). Another frequently used parameterization of this matrix element, motivated by the heavy quark expansion, reads where v = p/m Bs , v = k/m Ds , and the relative velocity w = v ·v = (m 2 Bs +m 2 Ds −q 2 )/(2m Bs m Ds ).
From the comparison of eqs. (4) and (7) one gets The form factor G(w) used in experimental analyses of the B (s) → D (s) ν decay is proportional to f + (q 2 ), and reads where, for convenience, we introduced Our main target is the determination of G(1), and therefore we are particularly interested in the dominant h + (1) term that can be easily obtained from f 0 (q 2 max ), Unfortunately, however, H(1) is not directly accessible from the lattice. Instead, we need to compute the form factors f 0,+ (q 2 ) at several small values of the D-meson three-momentum and then extrapolate H(w) to H(1). As we shall see in the following, we manage to work with w 1 but by staying very close to zero recoil and the uncertainty associated with this extrapolation is completely negligible. To be more specific, at | k| = 0, we compute so that where, for shortness, we write R 0 ≡ R 0 (q 2 (w)).

III. STRATEGY
To extract the form factors f 0,+ (q 2 ) from numerical simulations on the lattice, one first computes the correlation functions where the interpolating source operators, the pseudoscalar densities P cs and P bs , are sufficiently separated in the time direction so that for 0 t t S one can isolate the lowest lying states with J P = 0 − that couple to two source operators, and then extract the matrix element of the vector current between the two. The simplest choice would be the local operator, P hs =hγ 5 s, but for practical convenience one often resorts to the smearing technique that helps to significantly reduce the couplings to radially excited states. In other words, the lowest lying states are better isolated when smeared source operators are used and when the corresponding time interval, within which the matrix elements are extracted, becomes larger. As mentioned in the previous section, we also need to give the D s meson a few momenta | k| = 0, in order to study the behavior of H(w) as a function of w and extrapolate to H(1). To make those momenta small and remain close to the zero recoil point, it is convenient to compute the quark propagators, S q (x, 0; U ) ≡ q(x)q(0) , by imposing the twisted boundary conditions [15]. Those are easily implemented by rephasing the gauge field configurations according to, where U µ (x) stands for the gauge links, θ µ = (0, θ), and L is the size of the spatial side of the cubic lattice. The quark propagator computed on such a rephased configuration, can then be combined with the ordinary (untwisted) propagator, S q (x, 0; U ), into a two-point correlation function. The resulting lowest lying state extracted from the exponential fall-off has a three-momentum different from zero, | k| = | θ|π/L [15][16][17]. Then, by choosing θ = (1, 1, 1)×θ 0 , which also minimizes the discretization errors, one can tune θ 0 to an arbitrary small value and therefore explore the kinematical region of B s → D s ν decay very close to zero recoil, w 1 (q 2 q 2 max ).
For a given θ 0 , the continuum energy-momentum relation would give which is a good approximation for the small values of θ 0 chosen in this study. Otherwise one can use w = E(θ 0 )/m Ds , and the free boson dispersion relation on the lattice: In terms of quark propagators the correlation function (14) reads In addition to the above three-pont correlation functions, we also computed the two-point correlators that are necessary to remove the source operators and gain access to the vector current matrix element. From the large time behavior of the twopoint correlation functions we can extract m H and Z H = 0|hγ 5 s|H , where h (H) stands for either c (D s ) or b (B s ), and T is the size of the temporal extension of the lattice. With these ingredients we are able to extract the desired hadronic matrix element from decomposition and then by using the projectors (5) we get the form factors f 0,+ (q 2 ), as indicated in eq. (6). As already mentioned above, to make sure that the lowest lying states are well isolated, we employed the smearing procedure discussed in our previous works [18] where the full information concerning the smearing parameters can be found as well.
While working with the fully propagating b-quark on the lattice, i.e. without resorting to an effective theory approach, the most difficult problem is to deal with large discretization errors. The reason is mostly practical since the lattice spacings (a) accessible in current lattice QCD simulations are not small enough to satisfy m b a < 1. The charm quark, on the other hand, can be simulated directly and the discretization effects associated with its mass can be monitored by working at several small lattice spacings. Hence, the strategy is to perform the computations starting from the charm quark mass and then successively increase the heavy quark mass by a factor of λ so that after n + 1 steps one arrives at m b . For each value of the heavy quark mass m h = λ k+1 m c we compute the form factor G(1, m h , m c ), and evaluate the ratio of form factors computed at two successive heavy quark masses, while keeping other valence quarks and the lattice spacing fixed. In practice we compute where the first argument in the form factor is w = 1, the second is the heavy quark mass that we want to send to the physical b-quark, while m c and a are the charm quark mass and the lattice spacing, both of which are kept fixed. Each of these ratios can then be extrapolated to the continuum limit, lim The advantage of considering σ k (1) instead of G(1, λ k+1 m c , m c ) becomes apparent when considering the heavy quark mass dependence. In the continuum limit, thanks to the heavy quark symmetry, the form factor scales with the inverse heavy quark mass as where the non-perturbative coefficients g 0,1,2,... should be determined from the fit to lattice data. Keeping in mind the practical limitations that prevent us to ensure that the heavy quark masses are smaller than the inverse lattice spacing, it is clear that it is very challenging to disentangle the physical effects from lattice artifacts in g i 's, and in the dominant term g 0 in particular. Consequently the resulting G(1) ≡ G(1, m b , m c ) suffers from systematic uncertainty, the size of which is very difficult (if not impossible) to assess and therefore cannot be used for a precision determination of |V cb | from B(B (s) → D (s) ν). In contrast, the successive ratios of the form factors satisfy lim and therefore instead of extrapolating to the inverse b-quark mass, one actually interpolates to σ(m b ).
In the continuum limit, we then fit the lattice data determine s 1,2 , and interpolate to σ(m b ). Another interesting feature is that the expansion in inverse heavy quark mass is strictly valid in the heavy quark effective theory (HQET) and had we used eq. (22) we would have had to include the perturbative matching between our results (obtained in full QCD) to HQET, and then the result of extrapolation to the b quark mass should have been converted back to QCD. In the ratios of form factors, Σ k (σ k ), the matching to HQET and back becomes completely immaterial as the matching factors cancel to a large extent. We attempted including these corrections to our interpolation to σ(1, m b ), and the results remained would change by a few permil level only, thus completely immaterial for our purpose.
To get the physically relevant G(1), one starts from the elastic form factor, the value of which is by definition G(1, m c , m c ) = 1. The physically interesting B (s) → D (s) form factor is then obtained as a product of σ k (1) factors discussed above, namely In this study we choose n = 8, which gives where we used m MS c (2 GeV) = 1.14(4) GeV [19], and m MS b (2 GeV) = 4.91(15) GeV [14].

IV. LATTICE DETAILS
In this work we use the publicly available gauge field configurations that include N f = 2 dynamical light quarks, generated according to the twisted mass QCD action with maximal twist [21] by the European Twisted Mass Collaboration (ETMC) [20]. Main details about 13 ensembles of gauge field configurations are collected in tab. I. We computed all quark propagators by using stochastic sources, and then applied the so-called one-end trick to compute the correlation functions [20].
Since we use the smeared source operators, the factor |Z Ds | needed in eq. (20) depends on the momentum k given to the D s meson. We computed |Z Ds ( k)| for each of our θ 0 -values and after dividing out the correlation function C µ (t, k) in eq. (18), by the exponentials and couplings |Z Ds ( k)| and |Z Bs |, we looked for the plateau region to extract the desired matrix element [cf. eq. (20)]. After inspection, we fixed the plateaus to the intervals t β ∈ [10, 13] 3.8 , [10,13] 3.9 , [12,17] 4.05 , [17,19] in an obvious notation. Those plateaus are chosen to be common to all the sea quark masses considered at a given lattice spacing, and to all the heavy valence quark masses. Note that the threepoint correlation functions (14) are computed with t S = T /2 for all of our lattices. As an example, we illustrate in fig. 1 the quality of the plateaus corresponding to the matrix element D s ( k)|V i |B s ( 0) , and their sensitivity to the values of the threemomentum k used in this paper [or better, to the values of θ 0 in (18)]. In tab. I we gave the value of the charm quark mass in lattice units. Other heavy quark masses are simply obtained after successive multiplication by λ, with the physical µ b = λ 9 µ c . We note that the errors on the form factors become large for very heavy quarks. As mentioned in sec. III the computation of H(w) can be made only at w = 1. We tuned the values of the twisting angle θ 0 for each of our lattices in such a way as to make the corresponding w fixed [cf. eq. (17)]. More specifically, apart from the zero recoil point w = 1, we computed the form factors with four different non-zero momenta given to D s , corresponding to w ∈ {1.004, 1.016, 1.036, 1.062} .
Clearly, only a tiny extrapolation H(w) to H(1) is needed. Both a linear and a quadratic fit in w to reach H(1) lead to indistinguishable results, both results being small and further suppressed by the mass factor in eq. (9) so that G(1) is largely dominated by the form factor h + (1) evaluated according to eq. (11). For very heavy quarks (m h close to m b ) the effect of that extrapolation becomes visible, but since the form factor computed with such a heavy quark is dominated by discretization and large statistical errors they do not have any significant effect on our final results. We then computed the ratios of the form factors G(1) obtained at each two successive heavy valence quarks as indicated in eq. (21). Importantly, a strong cancellation of statistical errors leads to very accurate Σ k 's. The values of all Σ k 's are presented in tab. II. In fig. 2 we illustrate the situation for two values of k. From these plots we can see that our lattice data exhibit very little or no dependence on the light sea quark mass, nor on the lattice spacing. Note also that for larger heavy quark masses the errors on Σ k are larger, and therefore the corresponding continuum value σ k will have larger error as well.
The extrapolation of Σ k (1) to the continuum limit is performed by using the following form and then identify where m ud stands for the average of the physical up and down quark masses computed on the same lattices [19]. As anticipated from fig. 2 the values of β k and γ k , as obtained from the fit of our data to eq. (28), are consistent with zero. The resulting σ k (1) = σ(1, λ k+1 m c ) ≡ σ(1, m h ) are given in tab. III. Since our data do not exhibit a dependence on the sea quark mass we also attempted extrapolating Σ k (1) → σ k (1) by imposing β k = 0 in eq. (28). Results for the first few σ k (1) remain practically indistinguishable from those obtained by letting β k as a free fit parameter. For higher masses, namely for Σ 4−8 (1), the results of two continuum extrapolations remain compatible but the error bars in the case of a free β k are considerably larger. The problem of larger errors for which is illustrated in fig. 3. We then proceed as in eq. (24) and obtain The first (more accurate) result agrees with the only existing unquenched lattice QCD result, obtained for the light non-strange spectator quark [8].
To calculate our results in eq. (31) no renormalization constant was actually needed. This is convenient but not particularly beneficiary for our computation since the vector current renormalization constants have been computed non-perturbatively in ref. [23] to a very good accuracy (cf. values listed in tab. I). Therefore, we were able to perform several checks and instead of starting from the elastic form factor G(1, m c , m c ) = 1, we could have started from a k < n to compute G(1, λ k+1 m c , m c ) in the continuum limit, and then applied σ k+1 . . . σ n to reach the b-quark mass. For example, by using k = 3, in the case with β k = 0. This results is obviously completely consistent with the number given in eq. (31). To get the above result we also needed to perform a continuum extrapolation of G(1, λ 4 m c , m c ) by using the expression analogous to the one shown in eq. (28). We checked and observed that our lattice data for the form factor are also independent on the light sea quark when the valence quark masses are fixed, a behavior very similar to what is shown in fig. 2. Furthermore we checked that, after adding the cubic term in 1/m h to eq. (30), the resulting G(1) = 1.047(61), remains fully consistent with our main result given in eq. (31). Although th finite volume effects are not expected to affect the quantities computed in this paper, they could appear when the dynamical (sea) quark mass is lowered. In order to check for that effect we can compare our results obtained on the ensembles VI and VII which differ by the volume. The situation shown in fig. 2 is a generic illustration of the situation we see with all the other quantities: the form factors are completely insensitive to a change of the lattice volume. All these checks suggest that our result (31) obtained by using β k as a free parameter, remains stable and we take it for our final result, namely G(1) = 1.052(46).
Finally, we repeated the whole computation for the non-strange case, i.e. by keeping the sea and valence light quarks degenerate in mass. We obtained G(1) = 1.079 (29) and G(1) = 1.033(95), corresponding to β k = 0 and β k = 0 respectively. This latter number is not helpful in reducing the error bar of |V cb | extracted from B → Dµν decays. It shows, however, that the method employed in this work can be used to get a percent precision of G(1) even in the non-strange case provided the statistical quality of the data is substantially improved. Note also that our G(1) in eq. (33) agrees with the result obtained by the expansion around the BPS limit in ref. [24]. We end this discussion with a comment concerning the non-zero recoil situation (w = 1). The analysis is essentially the same as in the zero-recoil case described above. From the correlation functions (18) and by using the projector P + µ (5) we get the form factor f + (q 2 ) which is proportional to the desired G(w, λ k m c , m c ), cf. eqs. (8,9). The observations made in the analysis of G(1) concerning the independence on the light sea quark mass and on the lattice spacing remain true after switching from w = 1 to w = 1. The values are given in tab. IV, where we again report our results both in the case in which the parameter β k in the continuum extrapolation (28) is left free and in the case in which β k = 0 is imposed. The net effect in the latter case is that the resulting error is considerably smaller. Using the parameterization of ref. [5], which takes into account the relation between the curvature and the slope of G(w), namely , one could attempt to extract the slope ρ 2 from our data. Knowing that the window of w's we consider here is very short (27), a clean determination of ρ 2 would require very accurate values of G(w). In our case we only obtain ρ 2 = 1.2 (8), or in the case where we dismiss the dependence on the sea quark mass (when the errors on G(w) are smaller) we get ρ 2 = 1.1(3), both being consistent with the experimentally established ρ 2 = 1.19(4)(4) [25]. The same quality of result for ρ 2 is obtained if the data are fit to [26] G(w) G(1)

V. SCALAR AND TENSOR FORM FACTORS
Recently measured B(B → Dτ ν τ ) by the BaBar collaboration indicated about 2σ-discrepancy with respect to the SM estimate, obtained by combining the measured B(B → Dµν µ ) and the known information about the scalar form factor [27]. Since then, a number of studies appeared trying to explain that discrepancy by interpreting it as a potential signal of New Physics (NP) [28,29]. In the models with two Higgs doublets (2HDM), the charged Higgs boson can mediate the tree level processes, including B → D ν, and considerably enhance the coefficient multiplying the scalar form factor in the decay amplitude. For that reason it becomes important to get a lattice QCD estimate of f 0 (q 2 ). Furthermore, the model independent considerations of NP also allow a possibility of having a non-zero tensor coupling, in which case one more form factor appears in the decay amplitude. The tensor form factor f T (q 2 ) is defined via, where the renormalization scale dependence reflects the fact that the tensor density in QCD is a logarithmically divergent operator. In what follows the µ-dependence will be tacitly assumed. In this paper we report the result of the first lattice QCD computation of the tensor form factor in the region close to q 2 max . More specifically, in this section we compute which directly enter the expression for differential decay rate that can be found in eg. ref. [29]. Since the form factors f +,T (q 2 ) are not accessible at q 2 max (zero recoil, w = 1), we computed R 0,T (w(q 2 )) at w 1.
The extraction of R 0 (q 2 ) is practically straightforward. After applying the projectors P +,0 µ (5) to the matrix element extracted from the correlation functions (18), we combine them in the ratios R 0 (w, m h , m c ) = R 0 (w, λ k m c , m c ). Our goal is again to use the ratios of R 0 (w) computed at successive heavy quark masses, and then reach the point corresponding to the physically relevant R 0 (w, m b , m c ) through interpolation in inverse heavy quark mass. To this end, we first form where we indicate the momentum transfer w, the masses of quarks entering the weak vertex (m c , m h ), and the fact that Σ 0 k (w)'s are obtained at fixed lattice spacing, a. In tab. V we present the results for Σ 0 k (w) for a specific value of w = 1.016. Before discussing the heavy quark mass dependence we need to extrapolate to the continuum limit, lim a→0 Σ 0 k (w) = σ 0 k (w) by using a form similar to (28), thus assuming the linear dependence on the dynamical ("sea") quark mass and on the square of the lattice spacing. Since we work with maximally twisted QCD on the lattice, the leading discretization errors are proportional to a 2 [21]. After inspection, we found again that the dependence of the form factors on the sea quark mass is indiscernible from our data and that our results depend very mildly on the lattice spacing. Since the dependence on the sea quark mass is negligible, we again consider the continuum extrapolation by setting β k (w) = 0, separately from the case in which β k (w) are left as free parameters. The net effect is that the error on σ 0 k (w) = Σ k (w, 0, 0) is considerably smaller in the case with β k (w) = 0 and the data better respect the heavy quark mass dependence.
With several σ 0 k (w) in hands, we need to discuss the heavy quark interpolation. We first discuss its value in the infinitely heavy quark mass limit. Using the heavy quark effective theory mass formula m Bs,Ds = m b,c + Λ + (λ 1 + 3λ 2 )/m b,c [30], we can consider the ratio of form factors given in eq. (8), knowing that h + (w) scales as a constant with inverse heavy quark mass. One then deduces that, for the charm quark fixed to its physical value, Equivalently, R 0 (w, m h , m c ) = m h R 0 (w, m h , m c ), scales as a constant in the heavy quark mass limit, and the the corresponding σ 0 (w) is then be described by a form similar to eq. (23) and the physically relevant R 0 (w, m b ) could be obtained from We can also rewrite the above formula in terms of R 0 (w), as and therefore σ 0 k = σ 0 n /λ. In other words the interpolation formula to be used in this case is, where, again, our λ = 1.176. An illustration of that interpolation is provided in fig. 4 for one specific case of w. We see that our results obtained by assuming the independence of R 0 (w) on the sea quark scale better with the heavy quark mass than those obtained by letting the parameter β k (w) in eq. (39) free, although the two are compatible within the error bars.
It is interesting to note that these results are consistent with the values that can be obtained from the results quoted in recent literature (in the nonstrange case). More specifically, from the lattice results of ref. [12] one finds α = 0.020(1) GeV −2 , while from those reported in ref. [11] one finds α = 0.022(1) GeV −2 . Recent QCD sum rule analyses give α = 0.021(2) GeV −2 [31].
Note also that near zero recoil, the central value of our result R 0 (q 2 ) = 0.77(2), coincides with the quark model results of refs. [32,33].
To our knowledge, there is no QCD based determination of the B (s) → D (s) transition tensor form factor. The only existing result is the one presented in ref. [33] for the nonstrange case (B ud → D ud ν ) in which a constituent quark model has been employed. That obviously did not allow to keep track of the QCD anomalous dimension. However, as we shall see, their result [f T (q 2 )/f + (q 2 ) = 1.03(1)] is rather close to what we obtain from our lattice simulations. Furthermore, in ref. [33] it was found that this ratio R T (q 2 ) = f T (q 2 )/f + (q 2 ) is a flat function of q 2 .
On the lattice, the extraction of the form factor f T (q 2 ) is completely analogous to what we explained in the previous sections for f + (q 2 ) and f 0 (q 2 ). Heavy quark behavior of f T (q 2 ) is similar to that of f + (q 2 ), which is simple to see after applying the heavy quark equation of motion to the b-quark, 2 1+/ v 2 b = b, which in the b rest frame reads γ 0 b = b, and thereforecσ 0i b = −icγ i b, so that the heavy quark behavior of the form factor f T (q 2 ) resembles that of f + (q 2 ). We again define the ratios computed at two successive quark masses that differ by a factor of λ, which we then extrapolate to the continuum limit by using (45) Like in the previous cases, we observe that Σ T k (w) does not depend on the sea quark mass and its dependence on lattice spacing is insignificant within our error bars. For that reason we made the continuum extrapolation by imposing β k (w) = 0 and by leaving β k (w) as a free parameter. Results of that extrapolation, σ T k (w) ≡ σ T (w, m h , m c ) with m h = λ k+1 m c , are then interpolated in heavy quark mass to the b-quark, according to, which is shown in fig. 6. As in the case of R 0 (q 2 ), we also here need to extrapolate one of the ratios R T (q 2 ) to the continuum limit. We checked that for either k = 2, or 3 or 4 we end up with completely consistent results for The results given in tab. IV are obtained by choosing k = 3. Furthermore, we included the evolution of the tensor density from µ = 2 GeV, at which the renormalization constants have been computed, to µ = m b by using where, and where we used a s (µ) ≡ α s (µ)/π, for shortness. The anomalous dimension coefficients are known to three loops in perturbation theory and for N f = 4 in the MS scheme their values are [34]: where we also gave the first few β-function coefficients. Finally, in the computation we used Λ MS N f =4 = 296(10) MeV [35].
As it can be seen from tab. IV the error on f T (q 2 )/f + (q 2 ) is getting larger for larger values of w. We are therefore unable to check on the flatness of R T (q 2 ) that is valid in the infinitely heavy quark mass limit.

VI. SUMMARY AND PERSPECTIVES
In this paper we presented the results of our lattice QCD study of the exclusive semileptonic B s → D s ν decay form factors in the region near zero recoil (close to q 2 max ). The method employed here is the one proposed in ref. [13] that allows circumventing the problem of extrapolation in the inverse heavy quark mass and to reach the physical answer through interpolation. This is achieved by studying the successive ratios of form factors computed with heavy "b"-quark mass differing by a fixed factor of λ. In that way, in the continuum limit, these ratios have a fixed value for m b → ∞ and instead of extrapolating, one interpolates to the (inverse) b-quark mass.
We first computed the normalization to the vector form factor relevant to the extraction of the CKM matrix element |V cb | from B(B s → D s ν ) with the light lepton in the final state ∈ {e, µ}. We obtained and found that the method used here can also be employed to compute G(1) for the non-strange decay modes B(B → D ν ), provided the statistical sample of gauge field configurations were larger. We also observe that the above error bar can be significantly reduced if one imposed that the form factor ratios and the form factors themselves do not depend on the mass of the dynamical (sea) quark, which is essentially what we see with all of our lattice data (at all values of the lattice spacing). Thanks to the use of twisted boundary conditions imposed on the valence charm quark, we were able to explore the region of very small momenta given to D s , and therefore to compute the form factors for small recoil momenta w 1. Since we restrained our analysis to very small w's, we could not estimate the accurate value of the slope of the form factor G(w).
Instead, we computed the ratio of the scalar to vector form factors, R 0 (q 2 ) = f 0 (q 2 )/f + (q 2 ), that is needed to interpret the recent discrepancy between the experimentally measured B(B → Dτ ν τ )/B(B → Dµν µ ) and its theoretical prediction within the Standard Model. Since the scalar form factor contribution to the decay rate is helicity suppressed in the Standard Model, it is much more significant for the case of the τ -lepton in the final state than in the case of µ. This contribution is very important in various NP scenarios. In this paper we computed f 0 (q 2 )/f + (q 2 ) by using the same method of ratios and by restraining our attention to the small recoil region. Of several small w's, we quote Finally, in the models of physics beyond the Standard Model in which the tensor coupling to a vector boson is allowed, a third form factor might become important. Here we provide the first lattice QCD estimate of this (tensor) form factor f T (q 2 ) with respect to the vector one, f + (q 2 ). By employing the same methodology as above, in the MS renormalization scheme and at µ = m b we obtain We attempted repeating the same analysis for the case of the non-strange decay mode and found f 0 (q 2 0 )/f + (q 2 0 ) = 0.73(4) for q 2 0 = 11.6 GeV 2 (and f 0 (q 2 0 )/f + (q 2 0 ) = 0.75(2) if neglecting the dependence on the sea quark mass), in agreement with the results obtained in refs. [11,12]. Finally, in the non-strange case we get f T (q 2 0 )/f + (q 2 0 ) = 1.06 (12), which becomes 1.10(7) if neglecting the dependence on the sea quark mass. These values show that the prospects of using this method for computing the form factors for the non-strange decay modes, B → D ν , are promising provided the statistical quality of the data is improved.    Results of the continuum extrapolation of Σ k (1) to σ k (1) = σ(1, λ k+1 mc) ≡ σ(1, m h ) using eq. (28). Results of extrapolation with β k as a free parameter are shown separately from those in which the observed independence on the sea quark mass is imposed in the fit (28) by setting β k = 0.