B c → B s ð d Þ form factors from lattice QCD

We present results of the first lattice QCD calculations of B c → B s and B c → B d weak matrix elements. Form factors across the entire physical q 2 range are then extracted and extrapolated to the physical-continuum limit before combining with Cabbibo-Kobayashi-Maskawa (CKM) matrix elements to predict the semileptonic decay rates Γ ð B þ c → B 0 s ¯ l ν l Þ ¼ 52 . 4 ð 2 . 5 Þ × 10 9 s − 1 and Γ ð B þ c → B 0 ¯ l ν l Þ ¼ 3 . 10 ð 21 Þ × 10 9 s − 1 . The lattice QCD uncertainty is comparable to the CKM uncertainty here. Results are derived from correlation functions computed on MILC Collaboration gauge configurations with a range of lattice spacings including 2 þ 1 þ 1 flavors of dynamical sea quarks in the highly improved staggered quark (HISQ) formalism. HISQ is also used for the propagators of the valence light, strange, and charm quarks. Two different formalisms are employed for the bottom quark: nonrelativistic QCD (NRQCD) and heavy-HISQ. Checking agreement between these two approaches is an important test of our strategies for heavy quarks on the lattice. From chained fits of NRQCD and heavy-HISQ data, we obtain the differential decay rates d Γ =dq 2 as well as integrated values for comparison to future experimental results. DOI: 10.1103/PhysRevD.102.014513


I. INTRODUCTION
The semileptonic weak decays B þ c → B 0 sl ν l and B þ c → B 0l ν l proceed via tree-level flavor changing processes c → sW þ and c → dW þ parametrized by the Cabbibo-Kobayashi-Maskawa (CKM) matrix of the Standard Model. Associated weak matrix elements can be expressed in terms of form factors which capture the nonperturbative QCD physics. Precise determination of the normalization and the q 2 dependence of these form factors from lattice QCD will allow a novel comparison with future experiment to deduce the CKM parameters V cs and V cd . Lattice studies of other semileptonic meson decays that involve tree-level weak decays of a constituent charm quark include [1][2][3][4][5][6]. Precise determination of these CKM matrix elements is critical for examining the second row unitary constraint jV cd j 2 þ jV cs j 2 þ jV cb j 2 ¼ 1: This will complement other unitarity tests of the CKM matrix. It is possible LHCb could measure B þ c → B 0 sμ ν μ using Run 1 and 2 data. For example, normalizing by B þ c → J=ψμν μ would yield a constraint on the ratio V cs =V cb . Due to CKM suppression, a measurement of B þ c → B 0μ ν μ is likely to require many more B þ c decays. A lattice study of the B þ c → B 0 sl ν l and B þ c → B 0l ν l decays involves the practical complication of a heavy spectator quark. Care must be taken in placing such a particle on the lattice to avoid large discretization effects. We consider two formalisms for the b quark. A valence NRQCD [7,8] b quark, a formalism constructed from a nonrelativistic effective theory, is used to simulate with physically massive b quarks. A complementary calculation uses HPQCD's heavy-HISQ method [9][10][11]. Here, all flavors of quark are implemented with the HISQ [12] formalism. This is a fully relativistic approach which involves calculations for a set of quark masses on ensembles of lattices with a range of fine lattice spacings, enabling a fit from which the physical result at the b quark mass in the continuum can be determined. The method with an NRQCD bottom quark also uses HISQ for the charm, strange, and down flavors. This study will demonstrate the consistency of the NRQCD and heavy-HISQ approaches by comparing the form factors extrapolated to the physicalcontinuum limit.
In the limit of massless leptons, the differential decay rates for B þ c → B 0 sl ν l and B þ c → B 0l ν l are given by where V is the relevant associated CKM matrix element V cs or V cd and f þ is one of two form factors that parametrize the continuum weak matrix element hB sðdÞ ðp 2 ÞjV μ jB c ðp 1 Þi The 4-momentum transfer is q ¼ p 1 − p 2 , and only the vector part of the V − A weak current contributes since QCD conserves parity. The contribution of f 0 to the decay rate is suppressed by the lepton mass and hence irrelevant for the decays to eν e and μν μ . The phase space is sufficiently small to disallow decays to τν τ . Form factors are constructed from the matrix elements that are obtained by fitting the appropriate lattice QCD three-point correlator data. By calculating correlators at a range of transfer momenta on lattices with different spacings and quark masses, continuum form factors at physical quark masses are obtained and then appropriately integrated to offer a direct comparison with decay rates that could be measured in experiment.
In this study, we begin with Sec. II in which details of the lattice calculations are described. Section II A reports on the parameters and gauge configurations used to generate the propagators. Next, Sec. II B explains how the correlators are subsequently constructed for the two different treatments of the heavy spectator quark, as well as how the correlator data are fit to extract the matrix elements. Our nonperturbative renormalization method required to obtain the form factors is set out in Sec. II C. Section III presents results of the lattice calculations. Correlator fits are examined in Sec. III A, while Sec. III B discusses results for the renormalization of the local lattice vector current. In Sec. III C, the form factor data for the cases of an NRQCD spectator and a HISQ spectator are plotted alongside. Section IV is concerned with the methodology and results from fitting the form factor data. An extrapolation of the form factors to physical-continuum point is presented in Secs. IV D and IV E that shows how the form factors depend on the mass of the spectator quark. Finally, in Sec. V, we give our conclusions.

A. Parameters and setup
We use ensembles with 2 þ 1 þ 1 flavors of HISQ sea quark generated by the MILC Collaboration [13][14][15] and described in Table I. The Symanzik-improved gluon action used is that from [16], where the gluon action is improved perturbatively through Oðα s Þ including the effect of dynamical HISQ sea quarks. The lattice spacing is identified by comparing the physical value for the Wilson flow parameter ω 0 ¼ 0.1715ð9Þ fm [17] with lattice values for ω 0 =a from [18,19]. Our calculations feature physically massive strange quarks and equal mass up and down quarks, with a mass denoted by m l , with m l =m s ¼ 0.2 and also the physical value m l =m s ¼ 1=27.4 [5]. For sets 1-5 in Table I, strange propagators were reused from [20], a study of the pseudoscalar meson electromagnetic form factor. Light propagators were reused from [21], an extension of [20] to the pion. The valence quark masses used for the HISQ propagators on these gluon configurations are given in Table II. The valence strange and charm quark masses used here were tuned in [18,20], slightly away from the sea quark masses to yield results that more closely correspond to physical values. The propagators were calculated using the MILC code [22].
We work in the frame where the B þ c is at rest, and momentum is inserted into the strange or down valence quark through twisted boundary conditions [24,25] in the (111) direction. The values of the momenta used are given in Tables III and IV. The periodic boundary conditions of the fermion fields are modified by phases θ i , ψðn þ N xî Þ ¼ e iπθ i ψðnÞ ð 4Þ so that the usual lattice momenta q i ¼ 2πk i =aN x , for integers k i , are shifted by πθ i =aN x . The corresponding q 2 is then constructed by taking q 0 to be the difference in energies of the lowest lying initial and final states. The coefficients of operators corresponding to relativistic correction terms in the NRQCD action are given in Table III. The valence b quark masses used for the NRQCD propagators are also given there. The values were taken from [8], where the b quark mass was found by matching the experimental value for the spin-averaged kinetic mass of the ϒ and the η b to lattice data. For the calculation with an NRQCD spectator bottom quark, we use sets 1-5 in Table I. Bare heavy-quark masses am h used for the heavy-HISQ method are shown in Table IV. The selection of heavyquark masses follows [11]. As well as sets 3 and 5, the heavy-HISQ calculation makes use of a lattice finer than the five sets featuring in the calculation with an NRQCD spectator, set 6 in Table I. This is motivated by the necessity to avoid large discretization effects that grow with ðam h Þ (as ðam h Þ 4 at tree level) while gathering data at large masses that will reliably inform the limit m h → m b .

NRQCD spectator case
For the case of an NRQCD spectator quark, random wall source [26] HISQ propagators with the mass of the charm quark are calculated and combined with random wall source NRQCD b propagators to generate B þ c two-point correlator data. Two-point correlators for B 0 sðdÞ are generated similarly. The strategy of combining NRQCD random wall propagators and HISQ random wall propagators to yield two-point correlators was first developed in [27]. NRQCD propagators are generated by solving an initial value problem. This is computationally very fast compared to calculating rows of the inverse of the quark matrix.
The three-point correlator needed here is represented diagrammatically in Fig. 1. A HISQ charm quark propagator is generated by using the random wall bottom quark propagator as a sequential source. Following Appendix B in [28], and excluding a spacetime-dependent sign, the sequential source is given by the spin-trace where Γ is the gamma matrix structure at the operator insertion, S RW b is the random wall NRQCD propagator, and is the space-spin matrix which transforms the naive quark field to diagonalize the HISQ action in spin space.

HISQ spectator case
The case of a HISQ spectator quark proceeds similarly with the only difference being the use of a HISQ propagator instead of an NRQCD propagator for the bottom quark. Again, the charm propagator uses the spectator bottom quark propagator as a sequential source. Multiple masses are used for the spectator quark, each requiring a different charm propagator for the three-point correlator. Figure 2 shows the heavy-charm pseudoscalar meson masses that arise from calculations with the am h values in Table IV. On set 6, the finest lattice considered, we reach a value for M H c that is 80% of the physical B c mass.
The same strange and light random wall HISQ propagators on sets 3 and 5 are used in both the NRQCD and the heavy-HISQ calculations; thus, the data on these lattices in the two approaches will be correlated. However, the effect of these correlations is small in the physical-continuum limit since the heavy-HISQ data on sets 3 and 5 are the furthest away from the physical b quark mass point, and hence these correlations are safely ignored.

Fitting the correlators
The correlator fits minimize an augmented χ 2 function as described in [29][30][31]. The functional forms for the twopoint and three-point correlators follow from their spectral decomposition and include oscillatory contributions from the staggered quark time doubler. The matrix elements are related to the fit parameters V nn ½i; j through where J is the relevant operator that facilitates the c → sðdÞ flavor transition. The pseudoscalar mesons of interest are the lowest lying states consistent with their quark content, so we are only concerned with the matrix elements for i ¼ j ¼ 0 since we restrict E½k ≤ E½k þ 1 by using lognormal prior distributions for the energy differences. The presence of i, j > 0 terms is necessary to give a good fit and to allow for the full systematic uncertainty from the presence of excited states to be included in the extracted V nn ½0; 0. On each set, the two-point and three-point correlator data for both c → s and c → d at all momenta are fit simultaneously to account for all possible correlations. The matrix elements and energies are extracted and form factor values determined, along with the correlations between results at different momenta.

C. Extracting the form factors
The partially conserved vector current (PCVC) ward identity allows for a fully nonperturbative renormalization of the lattice vector current. Since the same HISQ action is used for the c and sðdÞ quarks that couple to the W þ in both the NRQCD and heavy-HISQ approaches, we have the PCVC identity relating the conserved (point-split) c → sðdÞ lattice vector current and the local lattice scalar density S. We choose a local lattice operator V μ loc ; thus, Eq. (9) must be adjusted by a single renormalization factor Z V associated with that operator, giving q μ hB sðdÞ jV μ loc jB c iZ V ¼ ðm c − m sðdÞ ÞhB sðdÞ jSjB c i: ð10Þ Since Z V is q 2 independent, in principle Z V need only be found at zero recoil where q μ has only a temporal component [3]. This avoids the need to calculate threepoint correlators associated with the spatial components of the vector current matrix element that appear in Eq. (10) for q ≠ 0. However, in practice, it is preferable to determine f þ near zero recoil through the spatial components of the vector current matrix element, albeit with the additional cost in computing three-point correlators with the corresponding insertion. As in [2], we combine Eqs. (3) and (10) to give a determination of f 0 solely in terms of the scalar density matrix element. We use Eq. (11) and calculation of the vector current matrix element to determine f þ and f 0 for the full q 2 range following [3,32]. Thus, we will calculate matrix elements of both the local scalar density J ¼ S and the local vector where V μ is the vector current matrix element, except at zero recoil where the denominator vanishes and f þ cannot be extracted. We find that using Eq. (12) near zero recoil is problematic since both the numerator and denominator grow from 0 as q 2 is decreased from the maximum value at zero recoil. For the case where the spectator is an NRQCD b quark, we instead use Eq.
This method gives much smaller errors near to zero recoil. Although mathematically equivalent to Eq. (12), extracting f þ through Eq. (13) does not suffer an inflation of error near zero recoil since both the numerator and denominator are nonzero for all physical q 2 . However, since V i appears explicitly in Eq. (13), three-point correlators with an insertion of V i need to be calculated. For the case of the spectator NRQCD b quark, the use of Eq. (13) is straightforward except that it requires inversions of the charm quark propagator from a different sequential source [see Eq. (5)] to allow for insertion of the current V i ¼ γ i ⊗ γ i in the mixed NRQCD-HISQ three-point function. Collecting V i at nonzero 3-momentum transfer in the NRQCD calculation will also test for any q 2 dependence of Z V that would appear as a discretization effect. Using Eqs. (11) and (12) or (13), form factor data at a variety of lattice spacings, light quark masses, and momenta are obtained from the energies and matrix elements.

NRQCD spectator case
For the case of an NRQCD spectator quark, the form factor extraction is complicated by the energy offset as a consequence of the subtraction of the b quark rest mass inherent in the NRQCD formalism. While physical energy differences are preserved with NRQCD quarks, energy sums are not. Consequently, Particle Data Group (PDG) [33] values are used where necessary. For example, we take when extracting the form factors. We use interpolating

HISQ spectator case
For the case of a HISQ spectator quark, we work only with local scalar and vector currents. Expressed in the spin-taste basis, we use γ 5 ⊗ γ 5 for the H sðdÞ interpolating operator and two different operators, γ 5 ⊗ γ 5 and γ 5 γ 0 ⊗ γ 5 γ 0 , for the H c interpolator. The first of these, γ 5 ⊗ γ 5 , makes a tasteless three-point correlation function when the scalar density operator 1 ⊗ 1 is used. The second, γ 5 γ 0 ⊗ γ 5 γ 0 , allows for a tasteless three-point correlation function when we use the local temporal vector current operator γ 0 ⊗ γ 0 [3]. This requires the calculation of two H c two-point functions with the two different choices of operator at both the source and the sink. The difference in masses between these two different tastes of H c meson is tiny and, although consistently taken care of, it has no impact on the calculation. This ratio is preferable to an effective energy defined using CðtÞ=Cðt þ 1Þ since the ratio in Eq. (15) better suppresses the oscillatory contributions in Eq. (7). Error bars are present in the figure but mostly too small to observe. We exclude t min =a data points from the beginning and end points of the correlators in our fits to reduce the contributions from excited states. For each of the cases of an NRQCD and HISQ spectator quark, we fit all of the correlator data to Eq. (7) on each set simultaneously to obtain the correlations between the fitted parameters. Consequently, the correlator fits involve a large covariance matrix. Without extremely large statistical samples of results, small eigenvalues of the covariance matrix are underestimated [34,35] and this causes problems when carrying out the inversion to find χ 2 . We overcome this by using an SVD (singular-value decomposition) cut; any eigenvalue of the covariance matrix smaller than some proportion c of the biggest eigenvalue λ max is replaced by cλ max . By carrying out this procedure, the covariance matrix becomes less singular. These eigenvalue replacements will only inflate our final errors; hence, this strategy is conservative. The SVD cut reduces the χ 2 =d:o:f reported by the fit because it lowers the contribution to χ 2 of the modes with eigenvalues below the SVD cut. In order to check the suitability of the SVD cut, we must test the goodness-of-fit from a fit where noise (SVD noise) is added to the data to reinstate the size of fluctuations expected from the modes below SVD cut, as described in Appendix D of [35]. The χ 2 SVD-noise =d:o:f is used to check the goodness of fit for both cases of spectator quark.

A. Correlators
Many fits were carried out with different SVD cuts, number of exponentials N, and positions t 2pt min and t 3pt min of the first time slice where the correlators are fit. We selected the fit of the correlators on each lattice for form factor extraction based on the χ 2 =d:o:f and Q value.
The parameters used in the fits of correlators with an NRQCD spectator quark are presented in Table V. The parameters given in bold are those used for our final fits. Other values are used in tests of the stability of our form factor fits to be discussed in Sec. IV B.
We fit the heavy-HISQ correlator data to Eq. (7) on each set simultaneously, including correlations between data with different values of twist, heavy-quark mass, and H s=d final state. Values for t min =a, the chosen SVD cut, the number of exponentials used in Eq. (7), and the resultant value of χ 2 =d:o:f including SVD noise are given in Table VI. We also include in Table VI fits

B. Vector current renormalization Z V
In this section, we give our results for the renormalization factor Z V for the vector current [Eq. (10)] and test for dependence of Z V on q 2 (for the case of an NRQCD spectator) and on the spectator quark mass (for the case of a HISQ spectator).
The vector current renormalization factor Z V computed at different momentum transfer with NRQCD b quarks shows no significant dependence on q 2 on each set, demonstrated by Fig. 5. Mild lattice spacing dependence is observed, however. For each momentum, we use the Z V found at the corresponding q 2 from Eq. (10).
The Z V factor in Eq. (10) is associated only with the local vector current operator and should be independent of the spectator quark. Z V values obtained in the different calculations are tabulated in Tables VII--IX. Good agreement is seen on set 5 at zero recoil between the results with NRQCD and heavy-HISQ spectator quarks. Dependence on the mass of the spectator quark is displayed in Fig. 6. The plot includes values from the analogous calculation for the D → K case [3]. For D → K, a charm quark decays into a strange quark, as in B c → B s , but here the spectator quark is a light quark, much less massive than the heavy spectator quark in B c → B s . The Z V from B c → B s and D → K in Fig. 6 are nevertheless in good agreement, demonstrating negligible dependence on the mass of the quark spectating the c → s transition.
It is also of interest to compare vector current renormalization factors for different masses of quark featuring in the current. For example, [36] calculates the localsγ μ s vector current renormalization factor from an η s → η s three-point correlation function at q 2 ¼ 0 on the 2þ1þ1 MILC ensembles. This gave very precise values and it was possible to fit Z V to a perturbative expansion in α s (including the known first-order term) along with discretization effects. The fit is plotted in Fig. 7 alongside Z V for c → s values determined in this study. This plot shows differing behavior as a function of a 2 . The value for Z V , determined nonperturbatively, is a combination of the underlying perturbative series in α s evaluated at a scale    (32) related to the lattice spacing and discretization effects that depend on how it was determined. Since the underlying perturbative series is common to different determinations, comparison will reveal the differing discretization effects. Figure 7 shows this in the comparison of our Z V values for the localsγ μ c current with those determined for the local sγ μ s current. In the limit of vanishing lattice spacing, where discretization effects vanish, the renormalization factors are in agreement. One might worry that the large errors appearing in Fig. 7 for thesγ μ c renormalization factors determined here would carry forward into our determination of the form factor f þ . However, the vector current matrix element at zero recoil, which contributes the dominant error in Z V , is highly correlated with the vector matrix elements at nonzero recoil. These correlations cancel in the ratio V 0 =V 0 ðq 2 max Þ appearing when using Eq. (10) to construct the renormalized current Z V V 0 appearing in Eqs. (12) and (13). Hence, the uncertainty in the renormalization factor is not a large contribution to our final uncertainty in the form factors. Figure 8 provides an example of the extracted values for the form factor f þ , comparing results from the NRQCD and heavy-HISQ spectator calculations. The lines on the figures connect data on the same set at a given am h value and are present as a guide only. The spread of the heavy-HISQ data for different heavy-quark masses is small, and the NRQCD and heavy-HISQ results are in good agreement on the fine lattice. Discretization effects are more noticeable for the case of an NRQCD spectator quark, especially on the coarsest lattices, sets 1 and 2. We believe that they result from the B c meson in the calculation since the effects are comparable to those seen in the B c meson decay constant study with NRQCD b quarks in [37]. Data points outside the physical region of momentum transfer are unphysical but nevertheless aid the fit.

A. z expansion
The four form factors, f 0 and f þ for each of the B c → B s and B c → B d processes, at all momenta on all the lattices, are fit simultaneously to a functional form which allows for dependence on the lattice spacing a and bare quark masses. FIG. 7. Z V from the s → s vector current from [36] (red crosses) and the c → s vector current from the heavy-HISQ calculation given here (green squares). The curve is the fitted perturbative expansion, including discretization effects, detailed in [36]. The red circle is an extrapolated value at the lattice spacing associated with the superfine lattice. . The NRQCD data are marked with circles, the heavy-HISQ data are marked with crosses, and finally the D → K values are given by diamonds. As expected, no significant dependence on the spectator mass is observed.
The fit is carried out using the lsqfit package [38] that implements a least-squares fitting procedure. As is now standard, we map the semileptonic region 0 < q 2 < ðM B c − M B sðdÞ Þ 2 to a region on the real axis within the unit circle through so that the form factors can be approximated by a truncated power series in z. Here, we choose the parameter t 0 to be 0 so that the points q 2 ¼ 0 and z ¼ 0 coincide. The parameter t þ is in principle the threshold for production of mesons, the lightest being D þ K, from the cs current in the t channel. It is convenient here, however, to work with t þ ¼ ðM B c þ M B sðdÞ Þ 2 , but this gives a very small range for z because then t þ ≫ t − . To correct for this, we rescale z.
The rescaling factor that we use is jzðM 2 p Þj −1 , where M p is the mass of the nearest cs or cd meson pole (we use the same mass for both vector and scalar form factors for convenience). For B c → B s , we take M p as the mass of the vector meson D Ã s and for B c → B d , the mass of D Ã0 . Thus, we define z p then has a range more commensurate to that for the corresponding D decay and the polynomial coefficients in z p are Oð1Þ. Coefficients of the conventional expansion in terms of z can easily be obtained from the expansion in z p . Using z p also avoids introducing large heavy mass dependence through the z transform in the heavy-HISQ case, which otherwise would require large Λ QCD =M H c coefficients in the heavy-HISQ fit. Note that in the case of the heavy-HISQ spectator, the B-meson masses above in t þ are replaced by the appropriate heavy meson masses at each value of am h (see Sec. IV C).

B. NRQCD form factor fits
The form factor results from the calculation with NRQCD spectator quark are fit to Here, the dominant pole structure is represented by a factor Pðq 2 Þ given by ð1 − q 2 =M 2 res Þ −1 with M res the mass of the relevant cs or cd meson (the vector meson for f þ and the scalar for f 0 ). We take the values of M res from current experiment [33]: GeV. We do not include uncertainties in these values since Pðq 2 Þ is a purely fixed factor designed to remove much of the q 2 dependence from the form factors. For our lattice results, uncertainties enter Pðq 2 Þ from the uncertainty in our determination of q 2 in physical units, including that from the determination of the lattice spacing. Pðq 2 Þ multiplies a polynomial in z p , and the polynomial coefficients are obtained from [39]. For the b quark, we take tuned values 1 of the quark mass from Table XII in [8].
For each of the sea and valence quark flavors, δm sea and δm val are given by giving estimates of the extent that the quark masses deviate from the ideal choices in which appropriate meson masses are exactly reproduced.
For prior values on the parameters in Eq. (19), we use 0(1) for A ðnÞ , B ðnÞ , and C ðnÞ , and 0.0(5) for κ ðjÞ . The power series in Eq. (18) is truncated to include up to the z 3 p term. Fits without a pole, i.e., Pðq 2 Þ ¼ 1, yield no statistically significant discrepancies. This is not surprising since the poles are far away from the physical region of q 2 , and so the pole effect on the form factor can be reasonably absorbed into the polynomial. Finally, the kinematic relation To ensure consistency, we convert values from [8] in lattice units to physical units by using the lattice spacing determined in [8] from the ϒð2S − 1SÞ splitting. is imposed on the fit as a constraint (we have tested that removing this constraint makes very little difference to the fit in fact and f þ ð0Þ − f 0 ð0Þ is zero to well within 1σ.).
Constraints on b ðnÞ from unitarity, as in the Bourrely-Caprini-Lellouch [40] and Boyd-Grinstein-Lebed [41] expansions, are unnecessary here since the full range of physical momentum transfer can be reached and so extrapolation in q 2 , which may benefit in accuracy from imposing these constraints, is not required. Hence, more complicated fit forms that impose additional physical constraints are not expected to be appreciably advantageous.
In Figs. 9 and 10, we demonstrate that the form factors in the physical-continuum limit are insensitive to the choice of the parameters in the fits of the correlators. As can be seen in the figures, the coefficients in the fits of the form factors are stable, within their uncertainties, as the correlator fits on different sets are varied.
The fitted form factors from the NRQCD spectator case exhibit errors no greater than 4% across the entire physical range of q 2 when tuned to the physical-continuum limit. Figures 11-14 show the results on all the lattices along with the fitted function for the form factors in the physicalcontinuum limit.
The z 0 p and z 1 p behavior of the form factors is well resolved by our fit to Eq. (19), as well as the ðam c =πÞ 2 z 0 p discretization effect. Table X summarizes the corresponding parameters from the fit. After fitting, other parameters FIG. 9. z p expansion coefficients, for the calculation with an NRQCD spectator quark, computed using the variations of correlator fit parameters listed in Table V Table V for the B c → B d form factors. The x coordinate is the same as that in Fig. 9.   FIG. 11. Lattice results and fitted f 0 form factor data for B þ c → B 0 sl ν l with an NRQCD b quark. The gray band shows the fitted form factor tuned to the limit of vanishing lattice spacing and physical quark masses.
show errors comparable to the width of their prior and are consistent with 0. In particular, quark mass mistuning coefficients simply return their prior value.

C. Heavy-HISQ form factor fits
We take a similar approach to fitting the form factor results for the case of a heavy-HISQ spectator. Now, we have results at multiple heavy-quark masses and the conversion from q 2 to z space [Eq. (16)] uses the values of M H c and M H s or M H d , as appropriate, from our calculation. We then rescale z at each m h as described in Sec. IVA [Eq. (17)]. This rescaling gives a similar z range for each m h and avoids introducing spurious dependence on m h that comes simply from the z transform.
The heavy-HISQ results are then fit to a form that is a product of Pðq 2 Þ and a polynomial in z p as for the NRQCD case. We now require a fit form for the polynomial coefficients that accounts for ðam h Þ 2n discretization effects as well as physical dependence on m h , however. Motivated by heavy quark effective theory (HQET), we express this physical heavy mass dependence as a power series in Λ QCD =M H c . The form factor data from the heavy-HISQ approach are fit to   where we only include the term proportional to δm val s for the B c → B s case. Pðq 2 Þ, δm, and the tuned masses have the same definitions as in the NRQCD case (Sec. IV B). In the physical-continuum limit, this form collapses to Pðq 2 Þ P n z n p A ðnÞ 000 . Again, we apply the constraint f 0 ð0Þ ¼ f þ ð0Þ in the continuum limit (by fixing A ð0Þ 000 to be the same in the two cases).
Results for the extrapolated form factors are given in Figs. 15-18 together with the corresponding lattice data. For the B c → B s case, a n and c n take prior values 0(1) and b n , d n , and e n take prior values 0(0.3) to reflect the fact that they enter through loop effects. In the B c → B d case, we take prior values of 0(1) for a n and e n and 0(0.3) for b n and d n . In both cases, we take prior values of 0(1) for A n ijk except for when i ¼ 1 or j ¼ 1 where we use a prior value of 0(0.3) to account for the HISQ one loop improvement.
As in the case of an NRQCD spectator quark, we present coefficients of the form factors fits from many different fits of the correlator data. Figures 19 and 20 show that the coefficients are insensitive to the choice of the parameters in the fits of the correlators.

D. Chained fit
The form factor functions tuned to the physicalcontinuum limit from NRQCD and heavy-HISQ are compared in Figs. 21-24 in z space. There is good agreement across the entire physical range of z, with particularly good agreement for the more accurate B c → B s case.
While the fit forms for the form factors from NRQCD and heavy-HISQ at Eqs. (18) and (23) differ in appearance, they both allow for effects of discretization and mistuning of the quark masses. In the continuum limit with physical masses, the two forms collapse such that the parameters A ðnÞ from Eq. (19) and A ðnÞ 000 from Eq. (23) coincide. Plotted among the functions from the heavy-HISQ and NRQCD calculations is a function arising from a "chained" fit where A ðnÞ 000 from the heavy-HISQ fit were used as prior distributions for A ðnÞ in the form factor fit forms in the NRQCD study. We label this fit NRQCD from heavy-HISQ in Figs. 21-24. As with the separate fits for each case of spectator quark, the form factors for B c → B s and B c → B d are fit simultaneously. This chained fit has χ 2 =d:o:f ¼ 1.3 and is consistent with both the separate fits. We make our final predictions for the decay rates and values for ΓjVj 2 using the chained fit. FIG. 19. z p expansion coefficients, for the calculation with a HISQ spectator quark, computed using the variations of correlator fit parameters listed in Table VI for Table VI for the B c → B d form factors. The x coordinate is the same as that in Fig. 19. We include the coefficients A ðnÞ 0;þ from the chained fit in the ancillary json file BcBsd_ff.json [42].

E. Dependence of the form factors on the spectator quark mass
In order to build up a picture of the behavior of form factors, it is interesting to ask: how do the form factors for c to s=d decay depend on the mass of the spectator quark? We can answer that question with our heavy-HISQ calculation because we have results at a range of spectator quark masses from m c upward (see Fig. 2). Our form factor fits (Sec. IV C) enable us to extrapolate up to m b . Our most accurate results are for the c to s decay case and we concentrate on that here. Figure 25 shows the fit curve from the heavy-HISQ results for f s þ and f s 0 as a function of the heavy-charm meson mass (as a proxy for the spectator quark mass). The form factor curves that are plotted are those for q 2 ¼ 0 (where f þ ¼ f 0 ) and for the zero-recoil point (q 2 max ). At q 2 max , the daughter meson is at rest in the rest frame of the H c meson. The q 2 value at q 2 max falls slowly as the heavyquark mass increases above m c because the mass difference between H c and H s mesons falls. Examining the region between M η c and M B c in Fig. 25 we see almost no dependence on the spectator mass. The form factor value that shows the most dependence is f þ ðq 2 max Þ. This is not surprising because f þ shows the biggest slope in q 2 close to q 2 max and hence sensitivity to the value of q 2 max . Note that the curve from the heavy-HISQ analysis agrees with the NRQCD results at a spectator mass equal to that of the b. As discussed in the previous subsection, the form factors obtained from the two calculations agree across the full q 2 range.
We can also investigate the behavior of the heavy-HISQ fit function as m h is taken below m c to m l where contact is made with results for D → K from [3]. For the form factors at q 2 ¼ 0, we have Pðq 2 Þ ¼ 1 and our fit form at Eq.
Here A, B, and ω n take prior values 0(2) and we do not include aΛ terms for data from [33]. We find this fit function reproduces our data, as well as the physical values, well. Figure 25 also shows the result of this downward extrapolation. While this extrapolation below m c is outside the region where HQET is expected to be valid, the curves nevertheless show approximately the correct amount of upward movement necessary to reproduce the D → K results in [3] for f þ and f 0 at zero recoil and q 2 ¼ 0. The form factors at q 2 ¼ 0 continue to show almost no spectator mass dependence, and this is in agreement with the D → K results.

F. Decay rate
The hadronic quantity required for determining the decay rate and branching fraction is the integral where V is the CKM element V cs or V cd .
where the CKM matrix elements are responsible for the first errors and the second errors arise from our lattice calculations. The dominant source of lattice QCD uncertainty is the fitting of two-point and three-point correlators described in Sec. II B 3. We can convert these results for the decay width into a branching fraction using the lifetime of the B c meson, 513.49(12.4) fs [43]. This gives where now the third uncertainty is from the lifetime.
We also present the ratio of the ΓjVj −2 for B c → B s to B c → B d taking correlations into account between the numerator and denominator. From the chained fit of B þ c → B 0 sl ν l and B þ c → B 0l ν l form factors, we obtain ΓðB þ c → B 0 sl ν l ÞjV cd j 2 ΓðB þ c → B 0l ν l ÞjV cs j 2 ¼ 0.809ð53Þ: In fact, the uncertainty is roughly the same as if we were to treat the numerator and denominator as uncorrelated.  and f þ (above) for B þ c → B 0 sl ν l in the physical-continuum limit, plotted against the entire range of physical q 2 . This fit is described in Sec. IV D.
FIG. 27. Final form factors from the chained fits of f 0 (below) and f þ (above) for B þ c → B 0l ν l in the physical-continuum limit, plotted against the entire range of physical q 2 . This fit is described in Sec. IV D.

V. CONCLUSIONS
We have reported here the first calculations of the decay rates ΓðB þ c → B 0 sl ν l Þ and ΓðB þ c → B 0l ν l Þ, demonstrating the success of lattice QCD in studying decays of heavylight mesons. The use of HISQ-HISQ c → sðdÞ currents allows for a nonperturbative renormalization using the PCVC. We used two different formulations for the spectator b quark, heavy-HISQ and NRQCD. Results from the heavy-HISQ calculations are in good agreement with the physical-continuum form factors derived from the calculations using NRQCD b quarks, giving us confidence in assessing and controlling the systematic errors in each formulation. Simulating at a variety of spectator masses in the heavy-HISQ calculation has provided a check of the spectator independence of the renormalization procedure for the vector current. The NRQCD study also accessed Z V away from zero recoil to scrutinize momentum independence.
Our final form factors from the chained fit that combines both NRQCD and heavy-HISQ results are plotted against q 2 in Figs. 26 and 27. The decay rates are predicted from our calculation with 4% and 6% uncertainty for ΓðB þ c → B 0 sl ν l Þ ¼ 52.4ð2.5Þ × 10 9 s −1 and ΓðB þ c → B 0l ν l Þ ¼ 3.10ð21Þ × 10 9 s −1 , respectively. There is scope for significant improvement should future experiment demand more precision from the lattice. Such improvement would be readily achieved by the inclusion of lattices with a finer lattice in the heavy-HISQ calculation. "Ultrafine" lattices with a ≈ 0.045 fm were used in [11] to provide results nearer to the physicalcontinuum limit with am h ≈ am b . Larger statistical samples could also be obtained on the lattices used here, at the cost of more computational resources.

ACKNOWLEDGMENTS
We are grateful to Mika Vesterinen for asking us about the form factors for these decays at the UK Flavour 2017 workshop at the IPPP, Durham. We are also grateful to Matthew Kenzie for discussions about the prospects of measurements by LHCb. We thank Jonna Koponen, Andrew Lytle, and Andre Zimermmane-Santos for making previously generated lattice propagators available for our use and Euan McLean for useful discussions on setting up the calculations. We thank the MILC Collaboration for making publicly available their gauge configurations and their code MILC-7.7.11 [22]. This work was performed using the Cambridge Service for Data Driven Discovery