Numerical approach to QED contribution in thelepton g − 2Masashi Hayakawa
Department of Physics, Nagoya University@Nara, 13-16 December, 2012

Purpose
I report the updated result for aµ (QED) (T.Aoyama, M.H T.Kinoshita
and M.Nio, Phys. Rev. Lett. 109, 111808 (2012))
1011 × aµ (QED) = 116 584 718.951 (80) .
The new ingredients are summarized as follows;(10)(1) compelete calculation of 5-loop (10th-order) coefficient, aµ(8)(2) improvement of precision of 4-loop (8th-order) calculation, aµ(3) recalculation of lower-order mass-dependent terms using thelatest values of mµ /me , mµ /mτ found in P. J. Mohr, B. N. Taylorand D. B. Newell, arXiv:1203.5425 [physics.atom-ph]

In particular, (1) plays the essential role in the substantial
reduction of its uncertainty;−11
δaµ (QED)[ignorance on a(10)µ ] = O(1) × 10

⇒ δaµ (QED) = 0.08 × 10−11 ≪ δaµ (next expr.) = O(1) × 10−11 .

Content• I introduce the basic features of the muon g − 2,
• I overview the current status of the muon g − 2, in particularof the QED contribution aµ (QED),• I explain a particular numerical approach to the high-order

perturbative QED calculation.
• The brief overview is given inT. Aoyama, M. H., T. Kinoshita and M. Nio,Prog. Theor. Exp. Phys. 2012, 01A107 (2012) .

What is g − 2 ?
• A charged massive particle ψ becomes magnetized in the formof magnetic dipole moment, if it has non-zero spinµψ = gψe~s,2mψ cwhere gψ is g-factor.
• If ψ has spin s = 21 , gψ |classically = 2.
• The deviation aψ ≡ (gψ −2) /2, called anomalous magnetic

dipole moment or “g − 2”, is a quantity predictable in anyrenormalizable quantum field theory.
• Here, we focus on the muon g − 2 (aµ ).
• The status for ae is described in backup slides.

Why is muon g − 2 ?
• The experimentally measured value, aµ (exp), of aµ consists of

the quantum-mechanical dynamics of the standard model,aµ (SM), and possibly those from extra structures, aµ (new);aµ (exp) = aµ (SM) + aµ (new) .
γStandard M odel
Extra Structuresµ• In order to explore the existence of aµ (new), we ask ifaµ (SM) differs from aµ (exp);aµ (exp) − aµ (SM) = aµ (new) 6= 0 ?

What is QED contribution ?

• Standard model contribution aψ (SM) (ψ = µ or e) consists ofthree types of contributions;aψ (SM) = aψ (QED) + aψ (QCD) + aψ (weak) .
• QED contribution, aψ (QED), is calculated by QED withphotons and charged leptons (e, µ, τ ) only.
• QCD contribution aψ (QCD) is calculated by QCD + QEDwith purely leptonic contribution (≡ aψ (QED)) subtracted.
• aψ (weak) consists of all the others, i.e., the diagrams with atleast one W boson, Z boson or Higgs boson.

Current status of muon g − 2
We have 2.8σ discrepancy between the measured value (aµ (exp))
and the theory of the muon g − 2 (aµ (SM));
1011 × aµ (exp) = 116 592 089 (63)1110

1011 × aµ (SM) = 116 591 840 (59)× {aµ (exp) − aµ (SM)} =249 (87)aµ is more sensitive to the still-unknown particle(s)/interaction(s)than ae , by (mµ /me )2 ≃ 40000, but also to QCD;
1011 × aµ (QCD) =6 967 (59) .

The uncertainty of aµ (SM) is now saturated by that of aµ (QCD),
due to the significant reduction of the uncertainty of aµ (QED) thisyear: 1011 × δaµ (QED) = 0.080 !

Current status of muon g − 2• New experiments (Fermilab(E989), J-PARC(E34)) are plannedto achieve the precision 1011 × δaµ (next exp) = O(1);1011 × aµ (LO. had. v.p.) =1110× aµ (NLO. had. v.p.) =6 949.1 (42.8)−98.4 (0.8)

1110

× aµ (had. lbyl) =116 (40)1011 × aµ (QCD) =6 967 (59)1011 × aµ (SM) = 116 591 840 (59)1011 × {aµ (exp) − aµ (SM)} =1011 × δaµ (next exp) =249 (87)
O(1)• aµ (LO. had. v.p.) requires knowledge on QCD with precisionof O(0.1) % .
Its evaluation mostly relies on the high precision experiment ofσ(e+ e− (s) → hadrons).

Current status on muon g − 2• It is indispensable to computeaµ (had. lbyl) ∼ O(100) × 10−11 ∼ O(aµ (exp) − aµ (SM)) bymeans of lattice QCD simulation .

QCDµ−• Improvement of precision of aµ (QED) and aµ (exp) assumes

development in aµ (QCD).

Requirement for aµ (QED) from δaµ (next exp)
• What should we do to realize

δaµ (QED) . δaµ (next exp) = O(1) × 10−11 ?
• To what order 2n of perturbative expansion is necessary toknow,∞ α nXa(2n)aµ (QED) =,µπn=1(2n)where QED predicts aµ ?
• Perturbative aspects of QED in aµ (QED) is quite different
from those in ae (QED);“ α ”2“ α ”3α+ O(1) ×+ O(1) ×πππ“ α ”4“ α ”5+ O(1) ×+ O(1) ×+ ··· .
ππ”““ α ”32αα+ O(10) ×aµ (QED) = 0.5 × + O(1) ×πππ“ α ”4“ α ”5+ O(100) ×+ O(1, 000) ×+ ··· .
ππae (QED) = 0.5 ×

Requirement for aµ (QED) from δaµ (next exp)
• The rough estimate of orders of magnitude will show that the10th order may also be relevant;aµ(8) × α 4= O(100) × 10−11 ∼ aµ (exp) − aµ (SM) ,π α 5(10)= O(1) × 10−11 ∼ δaµ (next exp) .
aµ ×π• These are exactly the reasons why(8)• we improve the numerical precision of aµ•(10)
we try to compute aµ ,over about 8 years !

(⇒ O(0.01) % !!!) ,2012 Update of aµ (QED)
Table: aµ (QED) at each order 2n, scaled by 1011 (T. Aoyama, M. H.,T. Kinoshita and M. Nio, Phys. Rev. Lett. 109, 111808 (2012) )order 2n
24
6810
sum

using α(Rb)116 140 973.318 (77)413 217.6291 (90)30 141.902 48 (41)381.008 (19)5.0938 (70)
116 584 718.951 (80)using α(ae )116 140 973.213 (30)413 217.6284 (89)30 141.902 39 (40)381.008 (19)5.0938 (70)
116 584 718.846 (37)(10)

The complete calculation of aµeliminates the uncertainty∼ O(1) × 10−11 ∼ δaµ (next exp), which has been present unless it is done.
Now, the uncertainties in aµ (QED) come mostly from1. statistical uncertainty in the Monte Carlo integration of the 8th-order terms,
2. uncertainty in the fine structure constant α (2n = 2).

Perturbative features of aµ (QED) and rough estimate
I next(2n)• estimate rough orders of magnitude of aµfor 2n ≥ 6,

 α 2 α 3α+ O(10) ×+ O(1) ×πππ α 4 α 5+ O(100) ×+ O(1, 000) ×+ ··· .
ππaµ (QED) = 0.5 ×• show the validity of our numerical result for the 8th and 10thorder terms;
aµ(8) = 130.879 6 (63) ,a(10)= 753.29 (1.04) .
µPerturbative features of aµ (QED) and rough estimate(2n)• Since mµ > me while mµ < mτ , aψ

is dominated by(2n)A2 (mµ /me ) in the following decomposition:
(2n)
(2n)
aµ(2n) = A1 +A2(2n)• The term A1






mµ(2n) mµ(2n) mµ mµ+A2,,+A3memτme mτ

is a pure number, called mass-independent(2n)(2n)
universally contributes to all al , and isterm. A1calculated by QED with electron only.
mµ(2n)• The term A2
represents the contribution of all
meFeynman diagrams with at least 
one e-loop but with nomµ(2n)τ -loop. Similarly for A2.
mτmµ mµ(2n)• A3,represents the contribution of all Feynman
me mτ
diagrams with both e-loop(s) and τ -loop(s).

Perturbative features of aµ (QED) and rough estimate
At the sixth order (n = 3), aµ receives large contribution throughthe light-by-light scattering due to virtual e− e+ , aµ (lbyl6e ).

e−µ−Perturbative features of aµ (QED) and rough estimate• The sixth-order light-by-light scattering contributionaµ (lbyl6e ) is given as follows; α 3 mµ× H ln+c ,aµ (lbyl6e ) ≃πme
mµ2ln∼ 5.33 , H = π 2 = 6.57973627 . . . ,me3c ≃ −2H (numerically) ⇒ {. . .} ≃ H × (5.33 − 2) ≃ 15 .
(The leading logarithmic approximation, ≃ 35, is valid up to afactor.)
We thus have
 α 3 α 3≃a(lbyl)=O(10)×× a(6).
µµ6eππPerturbative features of aµ (QED) and rough estimate• π 2 arises in the following manner;• For g − 2, O(α3 ) interactions must consist of two Coulombicinteractions (γµ ) and one hyperfine interaction (σµν q ν ).
• These Coulombic interactions act statically, yieldingiπδ(k 0 ) × iπδ(q 0 ).
• The enhancement of aµ (lbyl6e ) is caused by the dynamics

forming a bound state of a µ− and a e+ , called muonium(A.S.Yelkhovsky, Sov.J.Phys.49, 654 (1989)).

e−µ−Perturbative features of aµ (QED) and rough estimate• At the 8th-order,e−e−µ−will be dominant and is estimated as(8)

A2

„mµme

«

(6)

≃ A2

≃ 180 .

„mµ; l-lme

«

×

2ln3

„mµme

«

−59ﬀ×3Perturbative features“ of” aµ (QED) and rough estimate• Our computed result A(8)2mµme

= 132.685 2 (60) is roughly identical with theestimation ∼ 180.
• At the 10th order, the contribution of the diagrams (Set VI(a)) obtained byinserting two second-order electron-loop-induced vacuum polarizations in thesix-order light-by-light scattering diagrams will be dominant :(10)

A2

„mµme

«

(6)

≃ A2

„mµ; l-lme«

×2ln3

„mµme

«

−59ﬀ2

×6

≃ 1000 ,which gives(10)

aµ×“ α ”5π`´≃ 1000 × 6.76 × 10−14´`≃ 6.8 × 10−11δaµ (exp) = 63 × 10−11 .

• Our result obtained by complete calculation is consistent with rough estimate
(10)

A2

„mµme

«

[Set VI(a)] = 629.141 (12) .
«„mµ(10)= 742.18 (87) .
A2me

Numerical Approach to QED contribution• We employ the parametric integral formulation, whose basicpart was described in P.Cvitanovic and T.Kinoshita, Phys. Rev. D 10,3978 (1974).

• It deals with the integral on the Feynman parameter space.
• It intends to subtract UV divergence at the numerical level.
• The numerical subtraction is possible only if divergences aresubtracted in a pointwise way by the terms that are expressedas the integrands on the same Feynman parameter space;Z
XUVbare[dz] fG (z)−fF (z) ,F

where z = {zi } are Feynman parameters, and the sum istaken over all normal forests F of G.

Numerical Approach to QED contribution• The integral must be constructed separately for the individualvertex diagrams, or for a set of vertex diagrams (, which sharesimilar UV-divergent structure) which are related via a
Ward-Takahashi (WT) identity to a single self-energy-like
diagram G;
+WT

⇐⇒

+

≡G.

• The individual integral in general contains infrared (IR)divergence. The numerical calculation can be done for thequantity which are free from both UV and IR divergences;
ZXX∆MG = [dz] fGbare (z) −fEIR (z) .
fFUV (z)−FE

Numerical Approach to QED contribution• The UV subtraction terms are not the ones required from theon-shell renormalization condition. We thus have to add theresidual renormalization terms {Rk } to get the contribution to(2n)alfrom the gauge-invariant subset S of diagrams(IR-divergence cancels among them)XX(2n)al [Set S] =∆MG +Rk .
G∈SkOur construction of fFUV (z) and fEIR (z) guarantees that everyRk is given in terms of finite quantities at lower-order.
• The integration has been carried out for ∆MG and theconstituents of Rk with help of adaptive Monte Carlointegration routine, vegas, on RIKEN supercomputer
systems, RSCC and RICC.

Strategy for calculation of 10-th order QED contribution(8)(l = e, µ) required about 20 years(although it is corrected later on).

• The calculation of al• The number (12678) of Feynman diagrams at the 10th orderis 14 times larger than that (891) at the 8th order.
• Writing the numerical program for a 10th-order Feynman

diagram correctly is much harder than for the 8th-order.
• Thus, the 10th-order calculation was considered to need500 ∼ 1000 years to complete.
• We have completed the 10th-order calculation in less than 10years.

Strategy for calculation of 10-th order QED contribution• The most difficult gauge-invariant subset is Set V consistingof quenched-type (q-type) diagrams:

1. The number of q-type diagram at the 10th order is 6354.
2. A q-type diagram has complicated structure of UV and IR
singularities (, and thus requires many subtraction terms).
3. There are tenth-order q-type diagrams having linear IR
subdivergence. (The corresponding subtraction terms can nolonger be constructed just in the power-counting scheme.)

Strategy for calculation of 10-th order QED contribution• If we write 6354 programs manually, mistakes will bescattered over those programs randomly.
• If we write a code generator which produces 6354 programs(to calculate ∆MG for G ∈ Set V), mistakes will be strongly
correlated, and it is only necessary to manage the code
generator itself.
• The crucial point to successfully implement the code

generator is the invention of a systematic scheme, whichenables to subtract linear and higher IR singularities (T.Aoyama,M.H., T.Kinoshita and M.Nio, Nucl. Phys. B 796 (2008) 184.)• Quadruple precision is needed• to realize subtraction of linear IR divergence,(8)• to achieve the precision O(0.01) % for ae .

Strategy for calculation of 10-th order QED contribution• A code generator for q-type diagrams was implemented forarbitrary order 2n of perturbation. We have tested its validityfor 2n = 4, 6 and 8.
(8)• We found the incorrectness of the previous result for A1(T.Aoyama, M.H., T.Kinoshita and M.Nio, Phys.Rev.Lett. 99, 110406 (2007));(8)

A1 (old) = −1.7203 ,(8)A1 (new) = −1.9106 ,which affects to al ashi(8)al A1 (old) = −5.0313 × 10−11 ,hi(8)al A1 (new) = −5.5620 × 10−11 .
This causes a significant change to ae in light of its
experimental precision, δae (exp) = 0.028 × 10−11 .

”Modern history” of α−1(ae )109 × α−1(ae) − 137.035999
1000
800600[2]
400[1]
200

[3]

[4]

0−200−400−600−800[1] Washington U. exp ⇒ Harvard U. exp(8)[2] A1 was corrected (2007)[3] Improved Harvard U. exp[4] Theory updation (2012)

Strategy for calculation of 10-th order QED contribution• We have developed code-generators for several subsets ofdiagrams other than quenched-type diagrams.
• Some tenth-order diagrams were calculated analytically;• A. L. Kataev, Phys. Lett. B 284, 401 (1992) [Erratum-ibid. B710, 710 (2012)].
• S. Laporta, Phys. Lett. B 328, 522 (1994).
• J. -P. Aguilar, D. Greynat and E. De Rafael, Phys. Rev. D 77,093010 (2008).
• P. A. Baikov, K. G. Chetyrkin and C. Sturm, Nucl. Phys. Proc.
Suppl. 183, 8 (2008).

These works provided us to test the validity of our results forthe corresponding contributions.

Summary and future perspective• The calculation of 12, 672 number of Feynman diagrams at

the 10th order substantially reduced the uncertainty inaµ (QED), which is now well below the expected uncertainty
in the next-generation experiment.
• However, it is a reasonable question: Is our result correct ?
• Even the eight-order term (∼ aµ (SM) − aµ (exp)), has been

computed only by us.
• We have seen the consistency with rough orders of estimate for(8)(10)

aµ and aµ .
(8)
• Check for aµ by third persons is an important subject.
• The Harvard group is now preparing the new measurement ofae . Accordingly, the reduction of the uncertainty of ae (QED)
is strongly requested. We need
• optimized quadruple precision arithmetics,• more sophisticated adaptive Monte Carlo integrator(c.f. R.Arthur and A.D.Kennedy, arXiv:1209.0650 [physics.comp-ph] ).

Standard model prediction of electron g − 2The brand new value of the electron g − 2 is (Aoyama, M. H., Kinoshita,
Nio, Phys. Rev. Lett. 109, 111807 (2012) for ae (QED))
1012 × ae (QED) = 1159 652 180.07 (6)8th (8)10th (77)α(Rb)1012 × ae (QCD) =1.68 (3)1210 × ae (weak) =0.0297 (5)1210 × ae (SM) = 1159 652 181.78 (6)8th (8)10th (77)α(Rb) (3)QCD
1012 × ae (exp) = 1159 652 180.73 (28)Here
• The numerals in a parenthesis denote the uncertainty in thefinal few digits.
• The above uses α(Rb), which was obtained by the recent

determination of h/mRb via optical lattice technique(R. Bouchendira, P. Clade, S. Guellati-Khelifa, F. Nez and F. Biraben,
Phys. Rev. Lett. 106, 080801 (2011)).

Physical implication of electron g − 2• Since ae (QED) occupies 99.999 999 85 % of ae , thecomparison of ae (SM) with ae (exp) provides us with a
precision test of QED, at present.
• The precision of ae (QED) can be systematically improved inperturbation theory of QED;ae (QED) =∞Xn=1(2n)where QED predicts ae2n ≥ 4).

a(2n)×e α nπ,(with help of lepton mass ratios for8th-order QED contribution

I(a)

I(b)

I(c)

I(d)

II(a)

II(b)

III

IV(a)

IV(b)

IV(c)

IV(d)

V

II(c)

Figure: Typical vertex diagrams representing 13 gauge-invariant subsets
contributing to the eighth-order lepton g −2 (891 diagrams in total).

10th-order QED contribution

The new ingredients in the update of ae (QED) in 2012 are
1. completion of computation of 12, 672 number of 5-loop(10)
Feynman diagrams to get ae ,(8)2. improvement of numerical precision of ae ,3. improvement of mass-dependent terms at the lower orders
using the latest values of me /mµ , me /mτ found in P. J. Mohr,
B. N. Taylor and D. B. Newell, arXiv:1203.5425 [physics.atom-ph] .

Physical implication of electron g − 2• Due to the above efforts, the uncertainty in ae (SM) is now

dominated by that of α(Rb);1012 · ae (SM) = 1159 652 181.78 (6)8th (8)10th (77)α(Rb) (3)QCD ,1012 · ae (exp) = 1159 652 180.73 (28) .
• We suppose that no extra contribution exists in ae (exp).
• We can get α(ae ) by solving ae (SM, α) = ae (exp) withunknown α.
• Check of compatibility of α(ae ) with the others, such asα(Rb), α(Cs), α(q-Hall), α(JC), . . . , provides us with across-sectional understanding on various physical phenomena.