form factors from lattice QCD
Abstract
We report the first lattice QCD calculation of the form factors for the standard model treelevel decay . In combination with future measurement, this calculation will provide an alternative exclusive semileptonic determination of . We compare our results with previous model calculations, make predictions for differential decay rates and branching fractions, and predict the ratio of differential branching fractions between and . We also present standard model predictions for differential decay rate forwardbackward asymmetries and polarization fractions and calculate potentially useful ratios of form factors with those of the fictitious decay. Our lattice simulations utilize nonrelativistic QCD and highly improved staggered light quarks on a subset of the MILC Collaboration’s asqtad gauge configurations, including two lattice spacings and a range of light quark masses.
pacs:
12.38.Gc, 13.20.He, 14.40.Nd, 14.40.DfHPQCD Collaboration
I Introduction
The decay occurs at treelevel in the standard model via the flavorchanging chargedcurrent transition, making it an alternative to in the determination of from exclusive semileptonic decays. The difference in these processes, a spectator strange quark in vs a spectator down quark in , is beneficial for lattice QCD simulations, because it improves the ratio of signal to noise. Though this process has not yet been observed, its measurement is planned at LHCb and is possible during an run at BelleII. This provides a prediction opportunity for lattice QCD.
In addition to the calculation of form factors for , we also calculate their ratios with form factors for the fictitious decay. Such ratios are essentially free of our largest systematic error, perturbative matching. In combination with a future calculation of using a highly improved staggered (HISQ) quark, these ratios would yield a nonperturbative evaluation of the matching factor for the current with nonrelativistic QCD (NRQCD) quark. This matching factor would be applicable to and simulations using NRQCD quarks.
To include correlations among the data for both decays, correlation function fits must include vast amounts of correlated data. To make such fits feasible, we have developed a new technique, called chaining, discussed in Appendix A. In addition, the use of marginalization techniques developed in Ref. Hornbostel:2011 () significantly reduces the time required for the fits.
The chiral, continuum, and kinematic extrapolations are performed simultaneously using the modified expansion Na:2010 (); Na:2011 () with the chiral logarithmic corrections fixed by the results of hard pion chiral perturbation theory (HPChPT) Bijnens:2010 (); Bijnens:2011 (). The factorization of chiral corrections and kinematics, as found at oneloop order by HPChPT, suggests the modified expansion is a natural choice for carrying out this simultaneous extrapolation. We refer to the combination of HPChPT chiral logarithmic corrections and the modified expansion as the HPChPT expansion.
Ii Form factors and matrix elements
The vector hadronic matrix element is parametrized by the scalar and vector form factors
(1)  
where and . At intermediate stages of the calculation we recast in terms of the more convenient form factors
(2) 
where . In the meson rest frame, the form factors are simply related to the temporal and spatial components of the hadronic vector matrix elements,
(3)  
(4) 
The scalar and vector form factors are related to by
(5)  
(6) 
where is the kaon threemomentum. This discussion generalizes in a straightforward way for the matrix elements.
Iii Simulation
Ensemble  

C1  2.647(3)  0.005/0.05  0.8678  1200  2  0.0070  0.0489  2.650  0.28356(15)  
C2  2.618(3)  0.01/0.05  0.8677  1200  2  0.0123  0.0492  2.688  0.28323(18)  
C3  2.644(3)  0.02/0.05  0.8688  600  2  0.0246  0.0491  2.650  0.27897(20)  
F1  3.699(3)  0.0062/0.031  0.8782  1200  4  0.00674  0.0337  1.832  0.25653(14)  
F2  3.712(4)  0.0124/0.031  0.8788  600  4  0.01350  0.0336  1.826  0.25558(28) 
Ensemble averages are performed with the MILC Collaboration’s asqtad gauge configurations Bazavov:2010 () listed in Table 1. Valence quarks in our simulation are nonrelativistic QCD (NRQCD) Lepage:1992 () quarks, tuned in Ref. Na:2012 (), and highly improved staggered (HISQ) Follana:2007 () light and quarks, the propagators for which were generated in Refs. Na:2010 (); Na:2011 (). Valence quark masses for each ensemble used in the simulations are collected in Table 1 and correspond to pion masses ranging from, approximately, 260 MeV to 500 MeV.
Heavylight meson bilinears are built from NRQCD and HISQ quarks (for details see Ref. Na:2012 ()) and lightlight kaon (and similarly for the ) bilinears are built from HISQ light and quarks (for details see Ref. Na:2010 ()). From these bilinears we build two and three point correlation function data
(7)  
(8)  
(9)  
where indices specify quark smearing. We generate data for both a local and Gaussian smeared quark, with smearing function introduced via the replacement in Eqs. (7) and (9). Three point and daughter meson two point correlation function data are generated at four daughter meson momenta, corresponding to . In three point data, these momenta are inserted at in Fig. 1. The sum over in Eqs. (8) and (9) is performed using random wall sources with U(1) phases , i.e. . In the three point correlator a meson source is inserted at timeslice , selected at random on each configuration to reduce autocorrelations. The current is inserted at timeslices such that and the daughter meson is annihilated at timeslice . Prior to performing the fits, all data are shifted to a common . This three point correlator setup is depicted in Fig. 1. Additional details regarding the two and three point correlation function generation can be found in Ref. Bouchard:2013a ().
The flavorchanging current is an effective lattice vector current corrected through . The lattice currents that contribute through this order are
(10)  
(11) 
Matrix elements of the continuum vector current are matched to those of the lattice vector current according to
(12) 
where
(13) 
The matching calculation is done to one loop using massless HISQ lattice perturbation theory Monahan:2013 (). In implementing the matching, we omit contributions. Ref. Gulez:2007 (), which used asqtad valence quarks, found contributions of this order to be negligible. In Ref. Bouchard:2013a (), which used HISQ valence quarks, these contributions to the temporal component of the vector current were studied and again were found to be negligible. We also omit relativistic matching corrections. These, and higher order, omitted contributions to the matching result in our leading systematic error. An estimate of this error, and its incorporation in our fit results, is discussed in the following section.
Iv Correlation function fits
Two and three point correlation function fit Ansätze, and the selection of priors, closely follows the methods of Ref. Bouchard:2013a (). Two point data are fit to
(14) 
where tildes denote oscillating state contributions and is the simulated energy. The physical ground state mass is related to the simulation ground state energy by
(15) 
where GeV Gregory:2011 () is adjusted from experiment to remove electromagnetic, annihilation, and charmed sea effects not present in our simulations, and is the spinaveraged energy of states calculated on the ensembles used in the simulation and listed in Table 1. The quark smearing is indicated by indices . Kaon and two point correlator data are fit to an expression of the form^{1}^{1}1The zero momentum has no oscillating state contributions due to mass degeneracy of its valence quarks.
(16) 
Results of two point fits satisfy the dispersion relation and are stable with respect to variations in and the range of timeslices included in the fits, as demonstrated for kaon two point data in Ref. Bouchard:2013a ().
Three point correlation function data are described by
(17) 
where the three point amplitudes , , , and are proportional to the hadronic matrix elements. The ground state hadronic matrix element is obtained from
(18) 
where the factor of accounts for numerical factors introduced in the simulation and associated with taste averaging and HISQ inversion. In the correlator fits we include data for several temporal separations between the mother and daughter mesons. On the coarse ensembles we include data for while for the fine ensembles we include data.
On each ensemble we perform a simultaneous fit to two and three point correlation function data for the and decays, at all simulated momenta, including both spatial and temporal currents, and for the temporal separations listed above. This ensures correlations among these data are accounted for in the analysis. However, fits to such large data sets produce unwieldy data covariance matrices and are typically not convergent, or require a prohibitively large number of iterations. This can be partially addressed by thinning the data, e.g. by the use of singular value decomposition (SVD) cuts, but this reduces the accuracy of the fits.
To address this problem we introduce a technique, which we refer to as chaining, to simplify fits to very large data sets. Consider a data set consisting of correlators, . Before the fit, all fit parameters are assigned priors. Chaining first fits then uses the best fit mean values and covariances to replace the corresponding priors in subsequent fits. The updated set of priors is then used in the fit to . In this and all subsequent fits, correlations are accounted for between the data being fit and those priors which are best fit results from previous fits — this is an important step as it prevents “double counting” data. After this second fit, the priors are again updated according to the best fit mean values and covariances. This process is repeated for all correlators. The collection of best fit mean values and covariances following the fit to are the final fit results. Chaining is described in greater detail in Appendix A.
We combine the use of Bayesian Lepage:2002 (), marginalized Hornbostel:2011 (), and chained fitting techniques. Our final fit results use marginalization with a total of states accounted for, of which are explicitly fit. We refer to such fits with the shorthand notation, . States accounted for but not explicitly fit are marginalized in that their contributions are subtracted from the data prior to the fit. This technique reduces significantly the time required to perform the fits. In Fig. 2 we show the stability of the fits under variations in the numbers of states explicitly included and the total number of states accounted for in the fit.
Ensemble  

C1  0.8244(23)  0.7081(27)  0.6383(30)  0.5938(41) 
C2  0.8427(25)  0.6927(35)  0.6036(49)  0.536(12) 
C3  0.8313(29)  0.6953(33)  0.6309(30)  0.5844(46) 
F1  0.8322(25)  0.6844(35)  0.5994(43)  0.5551(56) 
F2  0.8316(27)  0.6915(38)  0.6119(43)  0.5563(61) 
Ensemble  
C1  2.087(16)  1.657(14)  1.378(13)  
C2  1.880(12)  1.412(16)  1.142(33)  
C3  1.773(11)  1.4212(84)  1.184(10)  
F1  1.878(13)  1.385(12)  1.158(13)  
F2  1.834(14)  1.396(10)  1.163(14) 
The form factor results from the correlation function fits are tabulated in Table 2 and additional details are given in Appendix B.
The form factors obtained from these fits preserve correlations resulting from shared gauge field configurations and quark propagators used in data generation. The preservation of correlations is demonstrated in the top panel of Fig. 3 where, e.g., significant correlations among the form factor fit results are seen at common momenta and nonzero correlations among form factors for the two decays is suggested. The bottom panel of Fig. 3 shows the distribution over all ensembles of correlations among form factors for the two decays. Accounting for these correlations is useful in our determination of the ratio of form factors for the two decays. Fit results for , and the resulting form factor ratios, are presented in Appendix D.
V Chiral, continuum, and kinematic extrapolation
The results of HPChPT Bijnens:2010 (); Bijnens:2011 () suggest a factorization, to at least oneloop order, of the soft physics of logarithmic chiral corrections and the physics associated with kinematics in the form factors describing semileptonic decays of heavy mesons,
(19) 
The logarithmic chiral corrections, calculated in Ref. Bijnens:2011 () for several decays, are independent of . An unspecified function characterizes the kinematics.
To obtain results over the full kinematic range one must include lattice simulation data over a range of energies. However, for any relevant physical scale (e.g. , ), at nominal lattice momenta and there is no convergent expansion of the unknown function in powers of . This is an inherent limitation of characterizing the kinematics in terms of energy. The energy of the daughter meson is a poor variable with which to describe the kinematics.
In contrast, the expansion Boyd:1996 (); Arnesen:2005 (); Bourrely:2010 () provides a convergent, modelindependent characterization of the kinematics over the entire kinematically accessible range. Combining a expansion on each ensemble^{2}^{2}2This assumes the general arguments on which the expansion is based hold for heavier than physical quark masses and at finite lattice spacing. with the HPChPT inspired factorization of Eq. (19) allows a simultaneous chiral, continuum, and kinematic extrapolation of lattice data at arbitrary energies. Because the chiral logs are the same for and , linear combinations (i.e. and ) factorize in the same way and have the same chiral logs. Motivated by these observations, we construct a HPChPTmotivated modified expansion, which we call the “HPChPT expansion”, and fit the lattice data of Tables 2 and 9, with accompanying covariance matrix, to fit functions of the form
(20)  
where [logs] are the continuum HPChPT logs of Ref. Bijnens:2011 (), and generic analytic chiral and discretization effects are accounted for by . Resonances above but below the production threshold, i.e. those in the range , are accounted for via the Blaschke factor, . Though not observed, we allow for the possibility of a state in , with choice of mass guided by Ref. Gregory:2011 (). Our fit results are insensitive to the presence of this state. The factorization suggested by HPChPT may not hold at higher order Colangelo:2012 () so we allow chiral analytic terms, which help parametrize effects from omitted higher order chiral logs, to have energy dependence (i.e. to vary with ).
We note that Eq. (20) is the modified expansion introduced in Refs. Na:2010 (); Na:2011 (), with the coefficients of the chiral logarithmic corrections fixed by the results of HPChPT. In the chiral and continuum limits
(21) 
and Eq. (20) is equivalent to the BourrelyCapriniLellouch parametrization Bourrely:2010 () of the form factors.
Following Ref. Bourrely:2010 () we impose a constraint on from the expected scaling behavior of in the neighborhood of . The resulting fit function for is
(22) 
We write , , and , explicitly exposing the dependence on and . This is useful in explaining the implementation of a second kinematic constraint we impose on the form factors. At the kinematic endpoint , the continuum extrapolated form factors and are equal, i.e. . We impose this constraint by fixing the coefficient ,
(23) 
Imposing this constraint results in the fit function for :
(24) 
In the fit functions for and , Eqs. (22) and (24), and [logs] are given by,
(25)  
(26)  
with implicit indices in Eq. (25) specifying the scalar or vector form factor. We account for momentumindependent and momentumdependent discretization effects in . The values of that enter the fit are the values from the simulation and are, of course, small. Finite volume effects in the simulation are included via a shift in the pion log Bernard:2002 (). The infinite volume limit is taken by setting this shift to zero. Eq. (26) gives the HPChPT Bijnens:2011 () result for the chiral logarithmic correction to form factors. These expressions make use of the dimensionless quantities
(27)  
(28)  
(29) 
where . We determine and on each ensemble using correlator fit results for meson masses and simulation momenta. Light and heavy quark discretization effects are accommodated for by making the mild functions of the masses, accomplished by the replacements
where is chosen so that as varies over the coarse and fine ensembles .
Lastly, we account for uncertainty associated with the perturbative matching of Sec. III. With the matching coefficients calculated in Ref. Monahan:2013 (), we find contributions to be of the total contribution to . Of this the majority, , comes from the one loop correction and from the NRQCD matching via . For we find contributions at this order to be , with coming from the correction and from the NRQCD matching. The matching error results from omitted higher order corrections, the size of which we estimate from observed leading order effects, where we conservatively use the larger 4%. Following the arguments outlined in Ref. Bouchard:2013a () we estimate the matching error to be the same size as the observed contributions and take the matching error to be 4%. This is equivalent to taking the matching coefficient to be four times larger than the matching coefficient (13 times larger than ). This uncertainty is associated with the hadronic matrix elements and therefore, by Eqs. (3) and (4), with and . To correctly incorporate it in the results for and we convert our fit functions for into , multiply by , where is a coefficient representing the matching error with a prior central value of zero and width 0.04, then convert back to before performing the fit. Schematically, we modify the fit functions, defined in Eqs. (22) and (24), by
(31)  
(32)  
(33) 
then we use to fit the results of the correlation function fits of Sec. IV. Conversions between the form factors and are performed using Eqs. (5) and (6).
The results of a simultaneous fit to the data for and , in which the maximum order of [specified by in Eqs. (22) and (24)] is 3 and , are shown relative to the data in Fig. 4 for . Details of prior choices and fit results are given in Appendix C.
We test the stability of this fit to the following modifications of the fit Ansätze:

Truncate the expansion at .

Truncate the expansion at .

Truncate the expansion at .

Drop momentumdependent and momentumindependent discretization terms in Eq. (25).

Drop the dependent discretization terms in Eq. (LABEL:eqdisc).

Drop the lightquark massdependent discretization terms in Eq. (LABEL:eqdisc).

Add the following nexttonexttoleadingorder (NNLO) chiral analytic terms to as defined in Eq. (25):
(34) 
Drop the sea and valencequark mass difference term from Eq. (25).

Drop the strange quark mistuning term from Eq. (25).

Drop finite volume effects, i.e. set in Eq. (26).
The stability of the fit results to these modifications is shown in Fig. 5, where results are shown at the extrapolated point. This point is furthest from the data region where simulations are performed and therefore is particularly sensitive to changes in the fit function. In Fig. 5 our final fit result, as defined by Eqs. (22) and (24) with and by Eqs. (25)–(LABEL:eqdisc), is indicated by the dashed line and gray band.
Modifications 1, 2, and 3 vary the order of the truncation in and demonstrate that by fit results have stabilized and errors have saturated. We therefore conclude that the error of the fit adequately accounts for the systematic error due to truncating the expansion.
Momentumdependent and momentumindependent discretization effects proportional to are removed in modification 4. This results in a modest increase in and a negligible shift in the fit result. This suggests our final fit, which includes the effects, adequately accounts for all discretization effects observed in the data.
In modifications 5 and 6 we remove heavy and lightquark massdependent discretization effects with essentially no impact on the fit. That our results are independent of lightquark mass dependent discretization effects suggests that staggered taste violating effects are accommodated for by a generic dependence.
Modification 7 tests the truncation of chiral analytic terms after nexttoleadingorder (NLO) by adding the NNLO terms listed in Eq. (34). This results in a slight decrease in but has no noticeable effect on the fit central value or error. From this we conclude that errors associated with omitted higher order chiral terms are negligible.
Differences in sea and valence quark masses, due in part to our use of HISQ valence and asqtad seaquarks, are neglected in modification 8. This results in a small increase in and negligible change in the fit results. We account for these small mass differences in our final fit, though this test suggests they are unimportant in the fit.
Effects due to strange quark mass mistuning on the ensembles are omitted in modification 9, resulting in a modest increase in and no change in the fit central value and error. We include these effects in our final fit.
Modification 10 results in nearly identical fit results, suggesting that finite volume effects are negligible in our data. We include these effects in our final fit results.
Vi Form Factor Results
In this section we present final results, with a complete error budget, for the form factors. We provide the needed information to reconstruct the form factors and compare our results with previous model calculations.
Fig. 6 shows the results of the chiral, continuum, and kinematic extrapolation of Sec. V, plotted over the entire kinematic range of . The form factors, extrapolated to , have the value .
vi.1 Fit errors for the HPChPT expansion
The inputs in our chiral, continuum, and kinematic extrapolation fits are data (the correlator fit results for and in Tables 2 and 9 with the accompanying covariance matrix) and priors. The total hessian error of the fit can be described in terms of contributions from these inputs, as described in detail in Appendix A. We group priors in a meaningful, though not unique, way and discuss the error associated with the chiral, continuum, and kinematic extrapolation based on these groupings. As the priors are, by construction, uncorrelated with one another, we can group them together in any way we find meaningful. The resulting error groupings are uncorrelated and add in quadrature to the total error. In Fig. 7 we plot the following relative error components as functions of :

experiment: This is the error in the fit due to uncertainty of experimentally determined, and other, input parameters. It is the sum in quadrature of the errors due to priors for the “Group I” fit parameters listed in Table 7. This error is independent of and subdominant.

kinematic: This error component is due to the priors for the coefficients in Eqs. (22) and (24). A comparison of the fit results from modifications 1, 2, and 3 in Fig. 5 shows that by the fit results have stabilized and errors have saturated. The kinematic error therefore includes the error associated with truncating the expansion. The extrapolation to values of for which we have no simulation data is controlled by the expansion. As a result, the growth in form factor errors away from the simulation region is due almost entirely to kinematic and statistical errors.

chiral: This error component is the sum in quadrature of errors associated with priors for in Eq. (25). These terms are responsible for extrapolating to the physical light quark mass and for accommodating for the slight strange quark mistuning and the small mismatch in sea and valence quark masses due to the mixed action used in the simulation. As shown in Fig. 7, these errors are subdominant and do not vary significantly with .

discretization: We account for momentumdependent discretization effects via the , and momentumindependent discretization effects via the , terms of Eq. (25). In addition we allow for heavy and lightquark massdependent discretization effects via the and terms in Eq. (LABEL:eqdisc). The discretization error component, which is essentially independent of , is the sum in quadrature of the error due to the priors for these fit parameters.

statistical: The statistical component of the error is due to uncertainty in the data, i.e. the errors from form factor fit results of Table 2. Simulation data exist for for and over the range for . Extrapolation beyond these regions leads to increasing errors.

matching: The matching error is due to the uncertainty associated with the priors for introduced in Eq. (32) and discussed in the surrounding text.
In addition to the largest sources of error, which we account for directly in the fit, there are remaining systematic uncertainties.
We simulate with degenerate light quarks and neglect electromagnetism. By adjusting the physical kaon mass () used in the chiral, continuum and kinematic extrapolation, we estimate the “kinematic” effects of omitting electromagnetic and isospin symmetry breaking in our simulation to be . It is more difficult to determine the size of the full effects. However, in general electromagnetic and isospin effects are expected to be subpercent. We assume the error in our form factor calculation due to these effects is negligible relative to other sources of uncertainty.
Our simulations include up, down, and strange sea quarks and we assume omitted charm sea quark effects are negligible. This has been the case for processes in which it has been possible and appropriate to perturbatively estimate effects of charm quarks in the sea Davies:2010 ().
Our final form factor results, multiplied by the Blaschke factor , are shown in Fig. 8 where they are compared with results from a model calculation using perturbative QCD (pQCD) Wang:2012 () and a relativistic quark model (RQM) Faustov:2013 (). Our results provide significant clarification on the form factors at large .
vi.2 Reconstructing Form Factors
Coefficient  Value  

0.315(129)  
0.945(1.305)  
2.391(4.671)  
0.3680(214)  
0.750(193)  
2.720(1.458)  
In the physical limit our form factor results are parametrized in a BCL Bourrely:2010 () form with coefficients [see Eq. (21)]. Including the kinematic constraint and terms through order , we have
(35)  
(36) 
where
(37)  
(38)  
(39)  
(40) 
and the resonance masses are and . The values of the coefficients , derived from the extrapolation fit results of Sec. V, and the associated covariance matrix, are given in Table 3. Note that it is necessary to take into account the correlations among the coefficients to correctly reproduce the form factor errors.
Vii Phenomenology
With the benefit of ab initio form factors from lattice QCD, we explore the standard model implications of our results. In this section we make standard model predictions for several observables related to the decay for and .
The standard model differential decay rate is related to the form factors by
(41) 
In Fig. 9 we plot predicted differential decay rates for and , divided by , over the full kinematic range of .