Calculation of the hadronic vacuum polarization contribution to the muon anomalous magnetic moment

We present a first-principles lattice QCD+QED calculation at physical pion mass of the leading-order hadronic vacuum polarization contribution to the muon anomalous magnetic moment. The total contribution of up, down, strange, and charm quarks including QED and strong isospin breaking effects is found to be $a_\mu^{\rm HVP~LO}=715.4(16.3)(9.2) \times 10^{-10}$, where the first error is statistical and the second is systematic. By supplementing lattice data for very short and long distances with experimental R-ratio data using the compilation of Ref. [1], we significantly improve the precision of our calculation and find $a_\mu^{\rm HVP~LO} = 692.5(1.4)(0.5)(0.7)(2.1) \times 10^{-10}$ with lattice statistical, lattice systematic, R-ratio statistical, and R-ratio systematic errors given separately. This is the currently most precise determination of the leading-order hadronic vacuum polarization contribution to the muon anomalous magnetic moment. In addition, we present the first lattice calculation of the light-quark QED correction at physical pion mass.

Introduction.-The anomalous magnetic moment of the muon a μ is defined as the deviation of the Landé factor g μ from Dirac's relativistic quantum mechanics result, a μ ¼ ½ðg μ − 2Þ=2. It is one of the most precisely determined quantities in particle physics and is currently known both experimentally (BNL E821) [1] and from a standard model theory calculation [2] to approximately 1=2 parts per million.
Interestingly, the standard model result a SM μ deviates from the experimental measurement a expt μ at the 3-4σ level, depending on which determination of the leading-order hadronic vacuum polarization a HVP LO μ is used. One finds [3][4][5][6]  , and a expt μ . This tension may hint at new physics beyond the standard model of particle physics such that a reduction of uncertainties in Eq. (1) is highly desirable. New experiments at Fermilab (E989) [7] and J-PARC (E34) [8] intend to decrease the experimental uncertainty by a factor of 4. First results of the E989 experiment may be available before the end of 2018 [9] such that a reduction in uncertainty of the a HVP LO μ contribution is of timely interest.
In the following, we perform a complete first-principles calculation of a HVP LO μ in lattice QCD þ QED at physical pion mass with nondegenerate up and down quark masses and present results for the up, down, strange, and charm quark contributions. Our lattice calculation of the lightquark QED correction to a HVP LO μ is the first such calculation performed at physical pion mass. In addition, we replace lattice data at very short and long distances by experimental e þ e − scattering data using the compilation of Ref. [10], which allows us to produce the currently most precise determination of a HVP LO μ .
Computational method.-The general setup of our nonperturbative lattice computation is described in Ref. [11]. We compute where fðq 2 Þ is a known analytic function [11] and Πðq 2 Þ is defined as P x e iqx hJ μ ðxÞJ ν ð0Þi ¼ ðδ μν q 2 − q μ q ν ÞΠðq 2 Þ with sum over space-time coordinate x and J μ ðxÞ ¼ i P f Q fΨf ðxÞγ μ Ψ f ðxÞ. The sum is over up, down, strange, and charm quark flavors with QED charges Q up;charm ¼ 2=3 and Q down;strange ¼ −1=3. For convenience we do not explicitly write the superscript HVP LO. We compute Πðq 2 Þ using the kernel function of Refs. [12,13] with CðtÞ ¼ 1 3 P ⃗ x P j¼0;1;2 hJ j ð⃗ x; tÞJ j ð0Þi. With appropriate definition of w t , we can therefore write The correlator CðtÞ is computed in lattice QCD þ QED with dynamical up, down, and strange quarks and nondegenerate up and down quark masses. We compute the missing contributions to a μ from bottom quarks and from charm sea quarks in perturbative QCD [14] by integrating the timelike region above 2 GeV and find them to be smaller than 0.3 × 10 −10 .
We tune the bare up, down, and strange quark masses m up , m down , and m strange such that the π 0 , π þ , K 0 , and K þ meson masses computed in our calculation agree with the respective experimental measurements [15]. The lattice spacing is determined by setting the Ω − mass to its experimental value. We perform the calculation as a perturbation around an isospin-symmetric lattice QCD computation [16,17] with two degenerate light quarks with mass m light and a heavy quark with mass m heavy tuned to produce a pion mass of 135.0 MeV and a kaon mass of 495.7 MeV [18]. The correlator is expanded in the finestructure constant α as well as Δm up;down ¼ m up;down − m light , and Δm strange ¼ m strange − m heavy . We write where C ð0Þ ðtÞ is obtained in the lattice QCD calculation at the isospin symmetric point and the expansion terms define the QED and strong isospin-breaking (SIB) corrections, respectively. We keep only the leading corrections in α and Δm f which is sufficient for the desired precision. We insert the photon-quark vertices perturbatively with photons coupled to local lattice vector currents multiplied by the renormalization factor Z V [18]. We use Z A ≈ Z V for the charm [19] and QED corrections. The SIB correction is computed by inserting scalar operators in the respective quark lines. The procedure used for effective masses in such a perturbative expansion is explained in Ref. [20]. We use the finite-volume QED L prescription [21] and remove the universal 1=L and 1=L 2 corrections to the masses [22] with spatial lattice size L. The effect of 1=L 3 corrections is small compared to our statistical uncertainties. We find Δm up ¼ −0.00050ð1Þ, Δm down ¼ 0.00050ð1Þ, and Δm strange ¼ −0.0002ð2Þ for the 48I lattice ensemble described in Ref. [18]. The shift of the Ω − mass due to the QED correction is significantly smaller than the lattice spacing uncertainty and its effect on CðtÞ is therefore not included separately. Figure 1 shows the quark-connected and quark-disconnected contributions to C ð0Þ . Similarly, Fig. 2 shows the relevant diagrams for the QED correction to the meson spectrum and the hadronic vacuum polarization. The external vertices are pseudoscalar operators for the former and vector operators for the latter. We refer to diagrams S and V as the QED-connected and to diagram F as the QEDdisconnected contribution. We note that only the parts of diagram F with additional gluons exchanged between the two quark loops contribute to a HVP LO μ as otherwise an internal cut through a single photon line is possible, which is part of a HVP NLO μ [23]. For this reason, we subtract the separate quantum-averages of quark loops in diagram F. In the current calculation, we neglect diagrams T, D1, D2, and D3. This approximation is estimated to yield an Oð10%Þ correction for isospin splittings [24] for which the neglected diagrams are both SU(3) and 1=N c suppressed. For the hadronic vacuum polarization the contribution of neglected diagrams is still 1=N c suppressed and we adopt a corresponding 30% uncertainty.
In Fig. 3, we show the SIB diagrams. In the calculation presented here, we only include diagram M. For the meson masses this corresponds to neglecting the sea quark mass correction, which we have previously [18] determined to be an Oð2%Þ and Oð14%Þ effect for the pions and kaons, respectively. This estimate is based on the analytic fits of (H7) and (H9) of Ref. [18] with ratios C m π;K 2 =C m π;K 1 given in Tab. XVII of the same reference. For the hadronic vacuum polarization the contribution of diagram R is negligible since Δm up ≈ −Δm down and diagram O is SU(3) and 1=N c suppressed. We therefore assign a corresponding 10% uncertainty to the SIB correction.  We also compute the OðαÞ correction to the vector current renormalization factor Z V used in C ð0Þ [18,20] and find a small correction of approximately 0.05% for the light quarks.
We perform the calculation of C ð0Þ on the 48I and 64I ensembles described in Ref. [18] for the up, down, and strange quark-connected contributions. For the charm contribution we also perform a global fit using additional ensembles described in Ref. [19]. The quark-disconnected contribution as well as QED and SIB corrections are computed only on ensemble 48I.
For the noisy light quark connected contribution, we employ a multistep approximation scheme with low-mode averaging [25] over the entire volume and two levels of approximations in a truncated deflated solver with all mode averaging (AMA) [26][27][28][29] of randomly positioned point sources. The low-mode space is generated using a new Lanczos method working on multiple grids [30]. Our improved statistical estimator for the quark disconnected diagrams is described in Ref. [31] and our strategy for the strange quark is published in Ref. [32]. For diagram F, we reuse point-source propagators generated in Ref. [33].
The correlator CðtÞ is related to the R-ratio data [12] by CðtÞ ¼ ½1=ð12π 2 Þ R ∞ 0 dð ffiffi ffi s p ÞRðsÞse − ffi ffi s p t with RðsÞ ¼ ½ð3sÞ=4πα 2 σðs; e þ e − → hadÞ. In Fig. 4 we compare a lattice and R-ratio evaluation of w t CðtÞ and note that the R-ratio data is most precise at very short and long distances, while the lattice data is most precise at intermediate distances. We are therefore led to also investigate a position-space "window method" [12,34] and write with ΔÞ − Θðt; t 1 ; ΔÞ, and a LD μ ¼ P t CðtÞ× w t Θðt; t 1 ; ΔÞ, where each contribution is accessible from both lattice and R-ratio data. We define Θðt; t 0 ; ΔÞ ¼ f1 þ tanh ½ðt − t 0 Þ=Δg=2 which we find to be helpful to control the effect of discretization errors by the smearing parameter Δ. We then take a SD μ and a LD μ from the R-ratio data and a W μ from the lattice. In this work we use Δ ¼ 0.15 fm, which we find to provide a sufficiently sharp transition without increasing discretization errors noticeably. This method takes the most precise regions of both data sets and therefore may be a promising alternative to the proposal of Ref. [35].
Analysis and results.- Table I shows our results for the window method and pure lattice determination. We quote statistical uncertainties for the lattice data labeled as (S) and the R-ratio data labeled as (RST) separately. Lattice and R-ratio uncertainties are added in quadrature. For the quark-connected up, down, and strange contributions, the computation is performed on two ensembles with inverse lattice spacing a −1 ¼ 1.730ð4Þ GeV (48I) as well as a −1 ¼ 2.359ð7Þ GeV (64I) and a continuum limit is taken. The discretization error (C) is estimated by taking the maximum of the squared measured Oða 2 Þ correction as well as a simple ðaΛÞ 4 estimate, where we take Λ ¼ 400 MeV. We find the results on the 48I ensemble to differ only a few percent from the continuum limit. This holds for the full lattice contribution as well as the window contributions considered in this work. For the quark-connected charm contribution, additional ensembles described in Ref. [19] are used and the maximum of the above, and a ðam c Þ 4 estimate is taken as discretization error. The remaining contributions are small and only computed on the 48I ensemble for which we take ðaΛÞ 2 as estimate of discretization errors.
For the up and down quark-connected and disconnected contributions, we correct finite-volume effects to leading order in finite-volume position-space chiral perturbation theory [36]. Note that in our previous publication of the quark-disconnected contribution [31], we added this finitevolume correction as an uncertainty but did not shift the central value. We take the largest ratio of p 6 to p 4 corrections of Tab. 1 of Ref. [37] as systematic error estimate of neglected finite-volume errors (V). For the SIB correction we also include the sizable difference of the corresponding finite and infinite-volume chiral perturbation theory calculation as finite-volume uncertainty. For the QED correction, we repeat the computation using an infinite-volume photon (QED ∞ [38]) and include the difference to the QED L result as a finite-volume error. Further details of the QED ∞ procedure are provided as Supplemental Material [39].
We furthermore propagate uncertainties of the lattice spacing (A) and the renormalization factors Z V (Z). For the quark-disconnected contribution we adopt the additional long-distance error discussed in Ref. [31] (L) and for the charm contribution we propagate uncertainties from the global fit procedure [19] (M). Systematic errors of the R-ratio computation are taken from Ref. [10] and quoted as Δm f ðtÞ to ðc 1 þ c 0 tÞe −Et , where we vary c 0 and c 1 for fixed energy E. The resulting p values are larger than 0.2 for all cases and we use this functional form to compute the respective contribution to a μ . For the QED correction, we vary the energy E between the lowest ππ and πγ energies and quote the difference as additional uncertainty (E). For the SIB correction, we take E to be the ππ ground-state energy.
For the light quark contribution of our pure lattice result we use a bounding method [42] similar to Ref. [43] and find that upper and lower bounds meet within errors at t ¼ 3.0 fm. We vary the ground-state energy that enters this method [44] between the free-field and interacting value [45]. For the 48I ensemble we find We quote the respective uncertainties as (E48) and (E64). The variation of ππ ground-state energy on the 48I ensemble also enters the SIB correction as described above. Figure 5 shows our results for the window method with t 0 ¼ 0.4 fm. While the partial lattice and R-ratio contributions change by several 100 × 10 −10 , the sum changes only at the level of quoted uncertainties. This provides a nontrivial consistency check between the lattice and the R-ratio data for length scales between 0.4 fm and 2.6 fm. We expand on this check in the Supplemental Material [39]. The uncertainty of the current analysis is minimal for t 1 ¼ 1 fm, which we take as our main result for the window method. For t 0 ¼ t 1 we reproduce the value of Ref. [10]. In Fig. 6, we show the t 1 dependence of individual lattice contributions and compare our results with previously published results in Fig. 7. Our combined TABLE I. Individual and summed contributions to a μ multiplied by 10 10 . The left column lists results for the window method with t 0 ¼ 0.4 fm and t 1 ¼ 1 fm. The right column shows results for the pure first-principles lattice calculation. The respective uncertainties are defined in the main text.   lattice and R-ratio result is more precise than the R-ratio computation by itself and reduces the tension to the other R-ratio results. Results for different window parameters t 0 and t 1 and a comparison of individual components with previously published results are provided as Supplemental Material [39].
Conclusion.-We have presented both a complete firstprinciples calculation of the leading-order hadronic vacuum polarization contribution to the muon anomalous magnetic moment from lattice QCD þ QED at physical pion mass as well as a combination with R-ratio data. For the former we find a HVP LO μ ¼ 715.4ð16.3Þð9.2Þ × 10 −10 , where the first error is statistical and the second is systematic. For the latter we find a HVP LO μ ¼ 692.5ð1.4Þð0.5Þð0.7Þð2.1Þ × 10 −10 with lattice statistical, lattice systematic, R-ratio statistical, and R-ratio systematic errors given separately. This is the currently most precise determination of a HVP LO μ corresponding to a 3.7σ tension a expt μ − a SM μ ¼ 27.4ð2.7Þð2.6Þð6.3Þ × 10 −10 ½4: ð7Þ The presented combination of lattice and R-ratio data also serves to provide additional nontrivial cross-checks between lattice and R-ratio data. The precision of this computation will be improved in future work by adding the missing subleading QED and SIB diagrams and simulations at smaller lattice spacings and larger volumes.
We would like to thank our RBC and UKQCD collaborators for helpful discussions and support. We would also like to thank Kim Maltman and Masashi Hayakawa for valuable discussions. We are indebted to Fred Jegerlehner  [2]. The innermost error bar corresponds to the statistical uncertainty.