Theory of the Anomalous Magnetic Moment of the Electron

: The anomalous magnetic moment of the electron a e measured in a Penning trap occupies a unique position among high precision measurements of physical constants in the sense that it can be compared directly with the theoretical calculation based on the renormalized quantum electrodynamics (QED) to high orders of perturbation expansion in the ﬁne structure constant α , with an effective parameter α / π . Both numerical and analytic evaluations of a e up to ( α / π ) 4 are ﬁrmly established. The coefﬁcient of ( α / π ) 5 has been obtained recently by an extensive numerical integration. The contributions of hadronic and weak interactions have also been estimated. The sum of all these terms leads to a e ( theory ) = 1 159 652 181.606 ( 11 )( 12 )( 229 ) × 10 − 12 , where the ﬁrst two uncertainties are from the tenth-order QED term and the hadronic term, respectively. The third and largest uncertainty comes from the current best value of the ﬁne-structure constant derived from the cesium recoil measurement: α − 1 ( Cs ) = 137.035 999 046 ( 27 ) . The discrepancy between a e ( theory ) and a e ( (experiment) ) is 2.4 σ . Assuming that the standard model is valid so that a e (theory) = a e (experiment) holds, we obtain α − 1 ( a e ) = 137.035 999 1496 ( 13 )( 14 )( 330 ) , which is nearly as accurate as α − 1 ( Cs ) . The uncertainties are from the tenth-order QED term, hadronic term, and the best measurement of a e , in this order.


Introduction
The anomalous magnetic moment of the electron was discovered in 1947 by P. Kusch and H. M. Foley [1].The Zeeman spectra of the gallium atom in a constant magnetic field were measured, and the gyromagnetic ratio (g value) of the electron was determined.If the orbital g value is assumed to be one as Dirac theory predicts, the spin g value of the electron is determined to be g = 2 × (1.001 19 ± 0.000 05). (1) The g value of the electron derived from the Dirac theory is exactly an integer two, and the difference between the measured g value and Dirac's two is called the anomalous magnetic moment of the electron: a e ≡ (g − 2)/2 = 0.001 19 (5). ( This tiny 0.1% deviation of the electron g, as well as the Lamb shift of the hydrogen atom discovered one year earlier [2] were the achievements of the state-of-the-art improvements in the frequency measurements at that time.Theory of charged particles and photons, named as quantum electrodynamics (QED), was also developed in the same era.It must be distinguished from its old version from the 1920s, which suffered from the divergence problem associated with perturbative calculation [3,4].The new idea, renormalization, enables us to compute finite physical quantities from the perturbation theory of QED with high precision.J. Schwinger showed that the value of the anomalous magnetic moment of the electron a e can be attributed to the one-loop effect of QED [5].Precision comparison between measurement and theory in 1947-1948 was the starting point of stringent tests of QED and the standard model of elementary particles.The very high precision achieved in testing a e is a consequence of small values of the fine structure constant α and the small electron mass.The Lamb shift could not serve as the high precision test of QED since it depends on other quantities such as hadronic interactions and hadronic masses.
In this article, we review the current status of the theory of the electron anomalous magnetic moment, including the contribution of Feynman diagrams of QED up to the tenth-order perturbation theory.Of course, our starting points are Feynman diagrams and corresponding integrals.However, these integrals, usually defined in the momentum space, are not convenient for numerical integration.Thus, we transform them to forms accessible to numerical integration.Currently, the numerical integration method is the only way to approach the tenth-order QED Feynman diagrams.The numerical computation methods of QED have been given in several references: our approach is given in [6][7][8][9][10][11][12], and alternatives are found in [13][14][15].In Section 2, the contributions to a e from all known sources are listed.In the subsequent sections, step-by-step instructions are given for the conversion of Feynman integrals to one accessible to numerical integration, working on the fourth-order diagrams as an example.

Summary of Contributions to a e
Before digging into the calculation method of a Feynman diagram, let us summarize the theoretical value of the electron anomalous magnetic moment a e .The standard model contribution to a e can be divided into three terms: a e = a e (QED) + a e (weak) + a e (hadron).
The QED contribution a e (QED) involves the leptons, which are electron (e), muon (µ), and tau (τ), and photons.The contribution involving the weak bosons (W ± , Z 0 ) is categorized into a e (weak).
The contribution from quarks or hadrons without a weak boson is denoted by a e (hadron).Since the electron is much lighter than the weak bosons and hadrons, the QED contribution dominates over others: a e (weak) and a e (hadron) amount to only 0.026 ppb [16][17][18][19] and 1.47 ppb [20][21][22], respectively, of the whole contribution.The QED contribution is further divided according to its lepton-mass dependence.Since the anomaly a e is dimensionless, lepton-mass dependence appears in the form of the ratio between lepton masses.Thus, we may rewrite a e (QED) = A 1 + A 2 (m e /m µ ) + A 2 (m e /m τ ) + A 3 (m e /m µ , m e /m τ ), (4) where m e , m µ , and m τ are masses of the electron, muon, and tau-lepton, respectively.The term A 1 is independent of the lepton-mass ratios and universal for all lepton species.Both mass-independent and mass-dependent terms can be calculated in a framework of the perturbation series expressed by Feynman diagrams.The expansion parameter of the perturbation series is the coupling constant of QED, namely the elementary charge e of leptons.Since it involves only even powers of e, the perturbation expansion is expressed in terms of the fine-structure constant α, which is proportional to e 2 .Thus, a Feynman diagram with n-loops appears in the 2n th -order of the perturbation series of QED.
In SI units, the fine-structure constant α is defined as where h is the Dirac constant, which is the Planck constant h divided by 2π, c is the speed of light, and 0 is the electric constant.As is well known, in the current SI system, c and 0 are the defined constants.
The new SI will be launched in May 2019, and the Planck constant h and the elementary charge e will become the defined constants, while the electric constant 0 will no longer be a defined constant.
In the natural units c = h = 0 = 1, the fine-structure constant has a much simpler expression: Hereafter, we shall use the natural units.Since the fine-structure constant α is a dimensionless constant, it has the same value regardless of the chosen units.
Since measurements find that α = 1/137.035 . .., we are lucky enough to have an effective perturbation series for QED calculation: Because of the renormalizability of QED, every A (2n) i is finite and calculable by using the Feynman-diagram techniques.The results of QED calculation are summarized in Table 1.By now, all coefficients A (2n) i up to the eighth-order have been found analytically [5,[23][24][25][26][27][28][29][30][31][32].They are consistent with numerical evaluations [33][34][35].The references quoted here are responsible only for the final stages of the calculation of each term, which were built upon much effort of many scientists over seventy years.Note that the mass-dependent terms of the sixth-and eighth-order terms, 2 , are not given in the closed forms.However, sufficiently higher order terms of the series expansion in a mass ratio were analytically obtained [30][31][32].Furthermore, the mass-independent eighth-order term A (8) 1 in [26] is not fully analytic.
One of the recent achievements in the QED calculation is the near-analytic result of the eighth-order mass-independent term A 1 obtained by S. Laporta [26].The Feynman integrals of the 891 eighth-order vertex diagrams are expressed in momentum space representation.They are decomposed into 334 master integrals (MI) according to the Laporta algorithm [36].These MIs are numerically evaluated with high precision, up to 9600 digits, using the difference equation and/or the differential equation methods.The numerical values are then fitted by using the PSLQ algorithm [37,38], which is a finder of integer relations, and the analytic expressions of the MIs are determined.The MIs for the eighth-order cannot be expressed only by the elementary functions.They contain harmonic polylogarithms and one-dimensional integrals of the products of elliptic integrals.The analytic expressions of several MIs related to the internal light-by-light scattering diagrams, shown as IV(d) of Figure 4, have not been determined yet.Thus, the result obtained by S. Laporta is near-analytic, but known up to 1100 digits.
The tenth-order terms have been calculated only numerically thus far [35,39,40], though some small sets of diagrams are known analytically.An independent numerical double check has not been carried out yet.The mass-dependent term of the tenth-order was reported in [39].
The mass-independent term, A 1 , has been calculated continuously on supercomputers even after some results were published in [40].The latest numerical result is consistent with 6.675 (192) given in Equation ( 16) of [40] and supersedes it.
Obviously, the mass-dependent terms A 2 and A 3 depend on the lepton-mass ratios provided by [41,42], and they are used as input parameters.However, they are the only input parameters of QED calculation.Notably, the mass-independent contribution A 1 is calculated without any input parameters.Principles, such as U(1) gauge theory, Lorentz and CPT symmetries, and renormalizability determine A 1 completely., where n denotes the 2n th -order of the perturbation theory of QED, are listed.No input parameter is used to compute A (2n) 1 , while the electron-to-muon mass ratio m e /m µ = 0.483 633 170 (11)×10 −2 [41] and the electron-to-tau mass ratio m e /m τ = 0.287 585 (19)×10 −3 [42] are used for The assigned uncertainties for the fourth-, sixth-, and eighth-order mass-dependent terms come from the lepton-mass ratios.The uncertainties of the tenth-order terms are due to numerical integration.The tau-lepton contributions to the tenth-order term have not yet been calculated, but they are suppressed by the factor (m µ /m τ ) 2 compared with the muon contributions.In order to obtain the theoretical prediction of a e , however, we need the input parameter determined from measurements of the nature.QED itself cannot determine what the fine-structure constant α is.Its value can only be derived from measurements.The quantum Hall resistance, which is named as the von Klitzing constant R K = h/e 2 , used to be the best method to determine the value of α.However, the current best is to use atomic recoil measurements.A quotient of the Planck constant and the mass of an atom X, h/m X , can be precisely determined by using atom interferometry [44].Two quantum states of an atom are spatially separated by transferring momenta through Bloch oscillation of the optical lattice made of laser beams [45].The obtained quotient is converted to

Value (Error) References
where R ∞ is the Rydberg constant and A r (X) and A r (e) are the relative atomic masses of an atom X and an electron, respectively, which is defined by m X /u (m e /u), u being the unified atomic mass unit.All three are precisely known and found in the CODATA 2014 adjustment [41].A r (X) is determined from the cyclotron frequency of an ion in the constant magnetic field.A r (e) is from the bound g factor of the electron, and R ∞ is from hydrogen spectroscopy.To determine the values A r (e) and R ∞ from the measured quantities, one needs the QED corrections that have been obtained from many theoretical calculations [41].In this sense, the value of α(X) determined from the quotient h/m X is not perfectly independent of QED.However, the uncertainties due to QED corrections are sufficiently small, and the uncertainty from the quotient h/m X dominates over others.Therefore, α(X) can be regarded as an independent determination of QED.
Recently, two new measurements on the hydrogen spectra have been carried out, and new values of R ∞ have been reported.One is the 2S-4P transition by A. Beyer et al. [46].Their R ∞ is different from the CODATA adjusted value by 3.3σ.Another is the 1S-3S transition by H. Fleurbaey et al., which reports a value of R ∞ consistent with CODATA [47].If we use the value of R ∞ from [46], it increases the value of α −1 by only 0.3 × 10 −8 , which is within the uncertainty of α due to the h/m X measurements.Thus, we use the CODATA2014 adjustment value for R ∞ .The values of α currently available are from the Rb atom [48] and Cs atom [49] measurements of h/m X and are obtained as α −1 (Cs) = 137.035999 046 (27).
Using them and hadronic and weak contributions to a e as listed in where the uncertainties are due to the fine-structure constant α, numerical evaluation of the tenth-order QED, and the hadronic contribution in this order.Note that the uncertainty of QED is now smaller than that of hadrons.The best measurement of a e was performed at Harvard.Their result was [50] a e (expt.)= 1 159 652 180.73 (28) The differences between theory and measurement are a e (theory : α(Rb)) − a e (expt.)=1.31 (77) a e (theory : α(Cs)) − a e (expt.)=0.88 (36) where the discrepancies are about 1.7σ and 2.4σ for α(Rb) and α(Cs), respectively.We can obtain one more value of α recalling that the theory of a e depends almost exclusively on α.Equating the theoretical formula of a e (3) with the measurement (14), we obtain α −1 (a e ) = 137.035999 1496 (13)(14)(330), (17) where the uncertainties are from the numerical evaluation of the tenth-order QED term, the hadronic contribution, and the measurement (14).The differences in α from the atomic recoil determinations are where the discrepancies are about 1.7σ and 2.4σ for α(Rb) and α(Cs), respectively.

QED Scattering Amplitude
In order to investigate the magnetic property of a single electron, let us consider the scattering of an electron by an external magnetic field.QED respects symmetries such as Lorentz, charge, parity, and time-reversal symmetries.Therefore, the scattering amplitude of an electron by an external photon field should be expressed only by two form factors F 1 and F 2 : where σ µν = i 2 (γ µ γ ν − γ ν γ µ ) and the external electron momenta are on-the-mass-shell: F 1 (q 2 ) and F 2 (q 2 ) are called charge (or Dirac) and magnetic (or Pauli) form factors, respectively.
If we follow the on-shell renormalization scheme so that the electron charge e is the very elementary charge observed in any non-relativistic systems, the charge form factor should be renormalized as On the other hand, there is no reason that F 2 (0) should disappear.In fact, a non-vanishing value F 2 (0) gives rise to an anomalous magnetic moment of the electron.Substituting the Gordon identity into (20), we can rewrite the scattering amplitude as In the non-relativistic limit and static limit of q → 0, only σ ij , i, j = 1, 2, 3 contributes, and the second term of ( 24) is reduced to where ψ is a two-component spinor and σ i is the Pauli matrix.The effective Hamiltonian of the interaction between a magnetic moment µ and an external static magnetic field B is where s = σ/2 is the spin of the electron.We can therefore identify the spin g-factor of the electron as Thus far, no successful explanation has been given for why the anomalous magnetic moment a e = F 2 (0) should be positive.We know that a e is positive because Schwinger's calculation of the one-loop Feynman diagram of QED turns out to be +1/2 as a coefficient of α/π.

Feynman Diagrams
In order to compute a e , we need to evaluate the scattering amplitude of the electron within the framework of QED perturbation theory.It is the scattering amplitude that we want to compute, so quantum corrections on external electron legs must be taken into account.However, if covariant gauge fixing for the photon field is chosen, self-energy corrections on external electron legs provide no effect after renormalization with the on-shell condition.Since we use the Feynman gauge for the photon field with the on-shell renormalization condition, it is sufficient for us to deal with one-particle irreducible (1PI) vertex diagrams.If a non-covariant gauge, such as the Coulomb gauge, is chosen, quantum corrections on external legs must be taken into account to compute a e [51].
The number of 1PI-vertex-Feynman diagrams of a given order of perturbation theory can be counted by using QED in zero-dimensional space and time [52,53].In zero dimensions, the photon and electron fields are just bosonic and fermionic numbers, respectively, and the path-integral over field variables can be exactly performed.Assume that QED consists of electrons and photons only.
The generating function Γ for the 1PI-vertex Green function is analytically calculable, and its closed form expressed by a modified Bessel function is given in [53].Its series expansion in the coupling constant e is found as The coefficient of the 2n th power of e is the number of 1PI-vertex diagrams of the 2n th -order perturbation theory of QED.Typical vertex Feynman diagrams of the second-, fourth-, sixth-, eighth-, and tenth-orders of the perturbation theory are shown in Figures 1-5, respectively.There are three main structures of these Feynman diagrams.One type is a diagram without an electron loop such as (e) in Figure 3, Group V in Figure 4, and Set V in Figure 5.The second type is a diagram with electron loops, but only of the vacuum-polarization (VP) type, such as (a)-(c) in Figure 3, Group I(a-d), II(a-c), III in Figure 4.The third type is a diagram including an electron loop to which four or more photons are externally attached.They, light-by-light scattering (LbyL) diagrams, appear for the first time in the sixth-order perturbation theory, as shown in (d) in Figure 3.The tenth-order diagram I(j) shown in Figure 5 contains a VP function consisting of two LbyL loops.This is classified as belonging to the VP type.
The computational difficulty of diagrams (of the same order) depends very much on their structures.Once the integral of a diagram without a fermion loop is constructed, insertion of VP functions is almost trivial, particularly when a VP is the second-or fourth-order VP function, which are known as compact and analytic forms.For the sixth-order and higher order VP functions, their construction based on Feynman diagrams is relatively easy.Since VP functions consist of closed loop diagrams, they are free from infrared (IR) divergence and suffer only from ultraviolet (UV) divergence.This is because IR divergence of a QED Feynman diagram occurs only when a massless photon is attached to on-shell electron lines.This never happens in a VP Feynman diagram.Numerical handling of UV divergences is much easier, especially for our recipe called the K-operation, than those of IR divergences.Furthermore, the Padé approximants for the sixth-and eighth-order vacuum-polarization functions consisting of a single fermion loop are available [54,55].They allow us to do an independent and rigorous check of the contributions involving the higher order VP functions.
Vertex diagrams including an LbyL loop start to appear in the sixth-order.This means that diagrams involving an LbyL loop have a simple structure of UV divergence.It is known that the Feynman integral of an LbyL diagram is particularly lengthy and complicated, and its analytic evaluation is very tough [6,26,56].Numerical calculation of them, however, does not face such problems at all.The simple divergence structure makes construction of the integral comparatively easy and convergence of numerical integration fast.
The most challenging is thus the computation of higher order diagrams without a fermion loop.In the remaining sections, we will discuss how to numerically calculate them.Other diagrams with fermion loops will be discussed elsewhere.Reprinted from [12].Reprinted from [12].

Ward-Takahashi Sum
Let us focus on QED vertex diagrams without a fermion loop.For the higher order terms of the perturbation theory, the number of Feynman diagrams grows factorially with the loop order [52].Furthermore, as the order increases, the complexity of the integral derived from a Feynman diagram drastically increases.Reduction of the number of diagrams, if possible, is thus highly desirable.To achieve this goal, we combine vertex diagrams sharing the same quantum corrections and relate them to a self-energy diagram through the Ward-Takahashi identity: where Λ µ (p, q) is the "sum" of vertex diagrams that are obtained by inserting an external photon vertex in a fermion line of the self-energy diagram Σ in every possible way.Differentiating both sides with respect to q ν and taking a vanishing momentum transfer limit q → 0, we obtain This equation enables us to obtain the magnetic moment amplitude of the 2n th -order diagram by projecting it out from the right-hand side of (30).Although they are not explicitly written, photon momenta are imposed on a cut-off Λ, and a small mass λ is given to a photon.These regularization parameters are safely removed after the finite amplitude is constructed.
The QED Feynman-Dyson rules introduced in many QED textbooks lead to a Feynman integral expressed in momentum space.Instead, we use the parametric Feynman-Dyson rules [7,9,10] that enable us to write down the Feynman integral in terms of Feynman parameters without referring to the momentum representation at all.The Feynman integral derived from a 2n th -order self-energy diagram can be expressed using Feynman parameters: z 1 , z 2 • • • and z a , z b , • • • are assigned from the left to the right to 2n − 1 electron and n photon propagators, respectively.All Feynman parameters are non-negative, and the sum of all is restricted to one.As an example, two self-energy-like diagrams of the fourth-order are shown in Figure 6.To illustrate the parametric Feynman-Dyson rules, we start with the Feynman-Dyson integral of the self-energy diagram M 4a of the left pane of Figure 6: where q 1 and q 2 are the momenta of photon propagators a and b, respectively, and p 1 = p + q 1 , p 2 = p + q 1 + q 2 , p 3 = p + q 2 are the momenta flowing on the electron propagators 1, 2, and 3.By introducing the Feynman parameters z 1 , z 2 , • • • and using the trick the denominators of all propagators are combined as The loop momenta q 1 and q 2 can be "diagonalized", and Equation ( 33) is rewritten as where The U function is a "Jacobian" from the momentum space to the Feynman parameter space.After the shift of the loop momenta, the numerator of the integrand of (31), F, is expressed as where k i 's are the linear combinations of the diagonalized loop momenta q d 1 and q d 2 and A i 's are functions of Feynman parameters.Because of its definition in (34), G is obviously related to A i as More systematic derivation of U and other functions of the Feynman parameters will be given in Section 4.3, which is based on graph theory and is easily applicable to any higher order diagrams.The right-hand side of (30) for any order of the perturbation theory can be obtained by performing several manipulations on the Feynman parametric representation of a self-energy diagram G.They are expressed as where (dz) G stands for integration variables over Feynman parameters subject to the constraint ∑ i z i = 1.The function U is a "Jacobian" for transformation from the momentum space to Feynman parametric space and V is the combined denominators of all propagators of the diagram G.
Explicit expressions of U and V will be determined later.The exact manipulations are given in [7,9,10].Briefly speaking, the procedure corresponding to the operator Z µν is to insert an external photon vertex to one of the electron propagators of a self-energy diagram and to differentiate this electron propagator with respect to the momentum q flowing in through the external photon vertex.
The procedure corresponding to the operator C µν is to insert an external photon vertex to one of the electron propagators and to differentiate with respect to q all other electron propagators to which an external photon vertex is not attached.The procedure corresponding to the operator N ν comes from the differentiation of the denominator function V with respect to p ν .The operator E ν is defined by E ν = ∂F/∂p ν , where F is the product of all numerators of the electron propagators of a self-energy diagram.
The procedures and their corresponding operators of γ-matrices are explicitly constructed for the self-energy diagram M 4a .There are three electron propagators in which an external photon vertex can be inserted.Thus, the Z µν -operator consists of three terms: where Z µν i is obtained by replacing the electron line factor l i ≡ / k i + A i / p + m of F in (36) by S µν i : C µν is constructed by picking up two electron line factors: where F µν ij for i < j can be obtained by replacing l i and l j of F by γ µ and γ ν , respectively.The coefficient C ij for i < j is defined in terms of B ij , which will be defined later in Section 4.3, as where Since the denominator function V contains the external momentum p as given in (34), the N ν -operator is trivially determined as where G is defined in (34), and its explicit form is given in (37).The explicit form of E ν -operator for M 4a is given by where F i is obtained by replacing l i of F by γ ν .

Projection of Anomaly Contributions
We are now ready to extract the anomalous magnetic moment contribution paying attention to its Lorentz structure.After projection operators are applied to (38) and (39), the magnetic moment amplitude has the form of where the functions E, C, N, Z are defined as and the projection operators of P 1 and P 2 are given by with m 2 = p 2 = 1.Hereafter, we set m 2 = p 2 = 1 and the photon mass λ = 0, unless we need to distinguish m 2 and p 2 .

Building Blocks B ij
We introduce the "correlation" functions B ij between two propagators i and j of the same diagram.The propagators i and j can be any of the electron or photon propagators.These B ij 's are the very basic building blocks of the Feynman parametric representation of a loop diagram.They are determined by and only by the topology of a loop diagram [57,58].
Let us introduce a chain diagram that is obtained when the external legs of a diagram are removed and no distinction is made between electron and photon propagators.Thus, chain diagrams of each order of the QED perturbation theory are identical to those of the scalar φ 3 theory.A line between two nodes of a chain diagram is associated with a sum of several Feynman parameters belonging to it, and its direction is freely assigned.For a given chain diagram with n loops with n > 1, the numbers of nodes and lines are 2(n − 1) and 3(n − 1), respectively.For n = 1, a chain diagram has no node, and it is just a circle.The number of topologically-distinguished chain diagrams at a given order of perturbation theory is quite small compared with those of Feynman diagrams.They are 1, 1, 2, 5, and 16 for the second-, fourth-, sixth-, eighth-, and tenth-order perturbation theory of QED.
A given n loop chain diagram has n independent closed circuits.In the case of a QED self-energy diagram without a fermion loop, independent closed circuits are easily identified.A closed circuit consists of a photon propagator and electron propagators that lie between the two vertices where this photon comes in and out.The direction of a closed circuit is chosen as the same as that of the electron propagators.
The diagrams M 4a and M 4b are reduced to the same chain diagram.For M 4a , the independent closed circuits c 1 and c 2 consist of two lines l 1 = z 1 + z a and l 2 = z 2 and two lines l 2 = z 2 and l 3 = z 3 + z b , respectively.For M 4b , the independent circuits c 1 and c 2 consist of two lines l 1 = z 1 + z 3 + z a and l 2 = z 2 and two lines l 2 = z 2 and l 3 = z b , respectively.
The n × 3(n − 1) loop matrix T iβ for a chain diagram is defined as follows: +1 if the line β is contained in the circuit c i and has the same direction of c i , −1 if the line β is contained in the circuit c i and has the opposite direction of c i , 0 if the line β is not contained in the circuit c i , where i indicates one of the independent chain circuits and β stands for a line number.The loop matrix and the same for M 4b .The n × n symmetric matrix U ij is derived from T iβ such that The Jacobian U in ( 38) and ( 39) is obtained as the determinant of the matrix U ij : The correlation function B αβ between two lines l α and l β of a chain diagram is defined by where U −1 kl is the inverse matrix of U ij .
The B ij between two propagators of a Feynman diagram is identical with B αβ , where the propagator i and j are contained in the lines α and β, respectively.The scalar current A i of an electron propagator introduced in (36) is expressed using B ij as The denominator function V is given by where For M 4a , we find U, B ij = B ji in terms of Feynman parameters as and for M 4b , It is easy to check that U of ( 58) and G of ( 37) with ( 55) and ( 59) are identical to those given in (35).
Even for a higher order diagram, the construction of the loop matrix T iβ is trivial.Once the loop matrix is formed, the building blocks of the Feynman parametric representation can be obtained.For instance, the loop matrix T iβ of the tenth-order diagram X 253 of Figure 7 is found as where the lines l β are chosen as and the closed circuits are numbered according to the lexicographical order of the contained photon propagators.With this matrix, Formulas ( 52)-( 57) and ( 43) lead to all necessary building blocks.

Integration over Loop Momenta
We perform integration over diagonalized loop momenta.In our formulation, this can be done by two steps, by which the loop momenta k i appearing in (36) should turn into "contractions".
Firstly, all k i 's in the numerator are replaced by the correlation functions.If a term of the numerator contains an odd number of k i 's, it is dismissed after integration.Up to the sixth-order diagrams, the numerator contains 0, 2, or 4 k i 's.Therefore, the contraction rules of two k i 's needed are where the coefficient 4 represents the space-time dimension.The extension to the eighth-or tenth-order diagrams is straightforward.Secondly, the powers of V are changed according to the number of contractions of k i 's.As is seen in ( 47), the inverse power of V starts with either n or n − 1.If the numerator originally contains 2m k i 's and is produced by m contractions of k i 's, the denominators change from V n or V n−1 to After applying the whole contraction procedures, the amplitude of the anomalous magnetic moment from the self-energy-like diagram G can be written in terms of the Feynman parameters as For n = 2, the unrenormalized amplitude M 4a thus has the form: and the explicit forms of the numerator functions are found to be For M 4b , the formal expression of the amplitude is the same as (67).Because A 1 = A 3 , B 12 = B 23 , B 11 = B 33 for M 4b , the numerator functions are much shorter and are given by (69)

UV Renormalization by K-Operation
Let us focus on the structure of UV divergence of an unrenormalized amplitude M (2n) G .Though the amplitude M (2n) G originates from a self-energy diagram, it does not suffer from an overall divergence associated with this 2n th -order self-energy diagram.This is because its UV divergence is dropped when the anomalous magnetic moment contribution is projected out from it.UV divergences of M (2n) G come from subdiagrams of a self-energy diagram G, either self-energy or vertex, of a lower order than the 2n th .
As is discussed in Section 3, the renormalization procedure we have to use is the on-shell renormalization condition.However, the on-shell renormalization is not suitable for numerical calculation of a Feynman integral.If we perform the on-shell renormalization for an unrenormalized amplitude, UV divergence can be canceled, but it brings new types of IR divergence into the integral.We, therefore, carry out the renormalization in two steps.
Our UV subtraction method is named as the K-operation, which is realized in practice as a simple power-counting rule.The form of a UV subtraction term is picked up by applying the K-operation to the unrenormalized amplitude of (66), which is expressed with the auxiliary functions U, V, A i , B ij , and C ij .The K-operation is also applied to these auxiliary functions.The resultant integral is written in the same Feynman parameter space of the unrenormalized amplitude, and point-wise UV subtraction is realized.With a little algebra on the Feynman parameters, the resultant integral for UV subtraction exactly decouples into the renormalization constant determined and the lower-order magnetic moment amplitude.This renormalization constant contains the same UV divergence of the on-shell renormalization constant, but it is free from IR divergence, unlike the on-shell one.
Let us examine the fourth-order cases.The diagram M 4a has a second-order vertex subdiagram S 12a consisting of the electron propagators 1 and 2 and the photon propagator a.To extract the UV divergent terms associated with this subdiagram, the K(1, 2)-operation is applied to M 4a .The K(1, 2)-operation in the Feynman parametric space is equivalent to taking the simultaneous limit of vanishing Feynman parameters of z 1 , z 2 , and z a , which form S 12a .
The K(1, 2)-operation on B ij and U reflects the decoupling of subdiagram S 12a and residual diagram G 3b ≡ G/S 12a as With these building blocks, A i , V, and G can be formally written in the same form as ( 55) and (56).
We also need to apply the K(1, 2)-operation to the amplitude expression (67).Among many terms of the numerator, only the terms explicitly proportional to B 12 can survive.It is obvious if we recall that B 12 is a Feynman parametric representation of the product of two loop momenta flowing in S 12a .From (67) and (68), the K(1, 2)-operated amplitude becomes After a little algebra on the integration variables of the Feynman parameters, the amplitude exactly decouples to lower order quantities as where L UV 2 is the term including the UV divergence of the second-order vertex renormalization constant L 2 .The numbers in the parenthesis of L UV 2 and M 2 in (72) indicate which electron propagators of the original diagram M 4a are involved in them.The L UV 2 is determined as the leading divergent component of L 2 , and thus the remainder L R 2 , as well, such that Note that L R 2 is completely free from UV divergence, but is IR divergent.M 2 in (72) is exactly the second-order anomalous magnetic moment and gives the Schwinger term α/(2π).
The diagram M 4b has a second-order self-energy diagram S 2b as a UV divergent subdiagram.The K(2)-operated building blocks are A i , V, and G are obtained in the same form of ( 55) and ( 56) with these building blocks.Since the subdiagram S 2b contains one electron propagator and no contraction occurs within S 2b , the K(2)-operated amplitude formally has the same form of M 4b : The numerator functions are obtained from (69) by letting z 2 = 0 and In addition to them, B 2i and A 2 in (69) must be replaced by A 2 B 1i and A 1 A 2 , respectively.The latter part of the K(2)-operation comes from the procedure of the wave-function renormalization.The numerator F of the self-energy diagram M 4b is given by The mass and wave-function renormalizations require that the term sandwiched between γ β and γ β of F is replaced by the mass and wave-function renormalization terms: where k 1 + A 1 p is the momentum flowing on the electron propagator adjacent to the self-energy diagram S 2b and δm 2 and B 2 are the mass and wave-function renormalization constants of the second-order.The replacement rules for the numerator functions by the K(2)-operation: are the consequence of multiplication of the adjacent momentum in front of the wave-function renormalization constant.For the K-operation of the higher order self-energy diagram, replacement rules are not so simple, but can be managed by changing the definition of V and A i slightly.Readers may consult [8][9][10].
After change of variables and integration by parts, we find where M 2 * is the second-order magnetic moment to which a two-point vertex is inserted and Note that B R 2 is free from UV divergence.

Forest Formula
To extend the K-operation procedure for the UV renormalization to a higher order diagram, we need to handle multiple divergences occurring from many subdiagrams contained in it.Let us first list all divergent subdiagrams in a given self-energy diagram.For instance, the tenth-order self-energy diagram X 253 shown in Figure 7 has five UV divergent subdiagrams: three subdiagrams of a self-energy type, S 3c , S 567de , and S 2345678bcde , and two subdiagrams of a vertex type, S 56d and S 67e .The K-operation renormalized amplitude of ∆ M X253 is formally given by where the electron line numbers of the K-operation are shown as its subscripts to compactify the notation.To construct the UV subtraction terms, the product (81) must be expanded, and the order of multiple K-operations must be changed according to the diagrammatic relations of subdiagrams.The relation of subdiagrams S α and S β can be classified into three cases: 1. disjoint: S α ∩ S β = φ S α and S β do not share an electron propagator.The operation K α K β exists, and K α and K β are commutable.

inclusion: S α ⊂ S β
All electron propagators of S α are also components of S β .The operation K α K β exists, and K α and K β are not commutable.K β must be applied first, and then, K α .
Expanding the product (81), we find that 23 UV subtraction terms are needed to make M X253 free from UV divergence [59]: Note that this is a simple realization of Zimmermann's forest formula [59].
For instance, we can write down the fourth-order amplitudes free from UV divergences: The term corresponding to K 12 K 23 is not present, because the relationship between S 12a and S 23b of M 4a is overlapping.The K-operation renormalized amplitude ∆ M 4a is finite.Since the K-operation acts on the integrand, not on the integral, pointwise UV subtraction is realized in the Feynman parameter space.Thus, ∆ M 4a is ready to go to numerical evaluation.We call it the finite amplitude and denote it as ∆M 4a without a prime.The amplitude ∆ M 4b is also free from UV divergence.However, it still suffers from IR divergence.We need to remove it using a pointwise subtraction method.

IR Divergence
There are two kinds of origins of IR divergence arising in the unrenormalized magnetic moment amplitude M G .Both are related to a self-energy subdiagram, but it is not the direct source of IR divergence.Suppose a self-energy-like diagram G has a self-energy subdiagram S. When the adjacent electron propagators of S become almost on-the-mass-shell, the residual diagram G/S yields IR divergence.Because we use the Ward-Takahashi sum for the magnetic moment amplitude, the self-energy subdiagram S may have two properties, either the self-mass δm S or the magnetic moment M S [11].
When the subdiagram S behaves as a self-mass δm S , the residual diagram G/S is a magnetic moment amplitude with an insertion of a two-point vertex M (G/S) * .This additional vertex increases the number of electron propagators and makes the IR behavior of the amplitude worse than the case without the insertion M G/S .
To avoid the IR divergence of this kind, we need to complete the mass renormalization with the on-shell condition.The UV divergent part of δm S has been already subtracted by the K S -operation.To subtract the remaining part of the self-mass contribution, we introduce a new procedure, the residual mass renormalization R such that where When the subdiagram S yields the magnetic moment M S , the IR behavior of the residual diagram G/S is similar to that of the vertex renormalization constant L G/S that is obtained by replacing S by an electron-photon vertex.Thus, we introduce one more procedure, the I-subtraction, such that where The UV divergences arising in the IR subtraction terms (85) and (87) are removed by using the K-operation renormalization.Decoupled products of (85) and ( 87) are merged into the same Feynman parametric space of M G by applying the technique to derive (72) or (79) inversely.
All IR divergences of the magnetic moment amplitude M (2n) G can be removed by the above two IR subtractions.For nested IR divergences appearing in a higher order diagram, we prepare the annotated forest formula and apply the R-and/or I-subtractions as needed.There are two types of relations between two self-energy subdiagrams S α and S β .Namely, One of the S's is a magnetic moment, and another is a self-mass: (b) Both are self-masses: R S α R S β .
(c) Two S's cannot simultaneously become magnetic moments: No double I-subtractions.(c) S α is a self-mass, and S β is a magnetic moment: R α I S β .
As an example, let us construct the IR subtraction terms for the tenth-order diagram X 253 shown in Figure 7.There are three self-energy subdiagrams S 3c , S 567de , and S 2345678bcde .Because δm R 2 = 0 in our definition, the residual mass renormalization for S 3c is not needed.We find The IR-subtraction terms for the unrenormalized amplitude M X253 are produced first by using the multiple R and/or I-subtractions.Then, the UV divergences of these IR-subtraction terms are removed by using the K-operation described in Sections 4.5 and 4.6.The latter part of the construction of the IR-subtraction terms is symbolically written as the application of R-and/or I-operations to the UV-free amplitude ∆ M X253 .
In a similar fashion, we apply the R-and I-subtractions to M 4b , and the finite amplitude ∆M 4b is obtained as This is ready to go to numerical calculation.In fact, because of our definition of the K-operation, we have δm 2 = δm UV 2 and δm R 2 = 0.This simplifies the calculation of the physical contribution of a e , especially for the higher order terms.

Residual Renormalization
The UV renormalization we have employed is not the standard on-shell renormalization.In addition, the K-operation renormalization, which is a simple power counting rule of contractions and Feynman parameters, violates the Ward-Takahashi identity between the renormalization constants: We also artificially added IR subtraction terms in order to make the amplitude numerically calculable on a computer.The residual and finite renormalization procedure is to be introduced to obtain the physical contribution to a e from the numerically-calculated finite amplitudes.During this process, the violation of the Ward-Takahashi identity of the renormalization constants is fixed.Then, the gauge invariance of the physical contribution to a e is guaranteed.It is also used to check the IR cancellation in the physical contribution to a e as stated by the Kinoshita-Lee-Nauenberg theorem [57,60].
Let us explicitly work out the fourth-order case.The standard renormalization procedure for the diagrams M 4a and M 4b can be expressed by using the finite amplitudes given in (83) and (90).We find Thus, the physical contribution from the gauge-invariant set of the fourth-order becomes where is guaranteed to be finite because of the Ward-Takahashi identity.Thus, all three terms in the right-hand-side of (94) are finite, and the physical contribution a 4 is free from IR divergences.After numerical integration, the finite quantities in (94) are found: and we obtain which is in agreement with the analytic result −0.344 166 • • • (α/π) 2 .
As you can see from the derivation of it, an explicit recipe of the K-operation is not essential to derive the finite formula (94).The on-shell mass renormalization is mandatory to make an integral free from IR divergence and available for numerical calculation.However, the separation of UV and residual terms of the vertex and wave-function renormalization constants can be arbitrary and is not necessarily the K-operation.For instance, S. Volkov used the separation such that the Ward-Takahashi identity holds [15,61]: where (i) indicates that the electron-photon vertex is inserted in the electron propagator i of the self-energy diagram S and the sum is taken over all electron propagators of S. Thus, no residual renormalization is required.In [15,61], this method is applied to numerical calculation of the eighth-order vertex diagrams without a fermion loop.The result is consistent with the previous numerical calculation formulated by using the K-operation [39] and also the analytic calculation [26].

Higher Order Calculation
For the higher order diagrams without a fermion loop, we have applied the same process described in Section 4 to calculate the contributions to a e .Because of the complexity and length of finite integrals of the higher orders, especially of the tenth-order, we automated many of the procedures described in Section 4 [10,11].
The automatic code generator of the integrand of a given self-energy-like diagram without a fermion loop is called GENCODEN.It is applicable to any order self-energy-like diagrams up to tenth-order and is extendable to even higher orders.From one-line information representing a self-energy diagram G of the 2n th -order, it generates the integrand of the finite amplitude ∆M G as a set of FORTRAN programs.It can be numerically evaluated with multi-dimensional integration algorithm, such as VEGAS [62].The workflow of GENCODEN is the following: 1.The one-line input of self-energy-like diagrams G is given.It is a sequence of the names of photons attached to the electron propagators from left to right.For example, the tenth-order diagram X 253 of Figure 7 is represented by "abccdedeba."2. The integrand of the unrenormalized magnetic amplitude M G is determined.More precisely, the numerator functions N i , E i , Z i ,, and C i are determined in terms of the building blocks and the scalar currents.3. The building blocks B ij 's and U are determined in terms of the Feynman parameters from the topology of G. 4. The forest formula for UV divergences of G is constructed.UV subtraction terms are then generated in terms of the building blocks.5.The UV limit of the building blocks B ij and U is taken for each UV subtraction term.6.The annotated forest formula for IR divergences of G is constructed.IR subtraction terms are then generated.The building blocks B ij 's and U of the UV subtraction terms that show the same decoupling of the subdiagrams are borrowed and used.7. The integrand of the finite amplitude ∆M G is constructed combining together all of the above.
The 6354 vertex diagrams of Set V of the tenth-order as shown in Figure 5 can be reduced to 389 self-energy-like subdiagrams.The 389 integrands, each consisting of about 100,000 lines of FORTRAN code, are generated by GENCODEN running on a personal computer.While numerical integration of the 389 finite amplitudes of ∆M G are being carried out on supercomputers, the residual renormalization formula of the 2n th -order is derived from the symbolic manipulation like (94).Neither numerical calculation, nor difficult algebraic calculation are needed at this stage.For the tenth-order case, we obtain The finite integrals ∆M 2n , ∆LB 2n , ∆δm 2n are obtained from the magnetic moment amplitudes, the sum of vertex and wave-function renormalization constants, and the mass renormalization constants, respectively, of the 2n th -order diagrams without a fermion loop.An asterisk ( * ) indicates that the quantity can be derived from diagrams having a two-point vertex insertion.For instance, ∆L 4 * is the sum of finite parts of 24 diagrams, which are fourth-order vertices with a two-point vertex in one of four electron propagators.Some quantities appearing in (99) are identical with those used in the lower order calculations such as the sixth-and eighth-order contributions.Agreement between the analytic results and the numerical results obtained by using these quantities is indirect, but strong evidence that they are correct.
Newly-evaluated ones specifically for the tenth-order are ∆LB 8 , ∆δm 6 , and ∆L 4 * .The sixthand fourth-order quantities are easily calculated.For ∆LB 8 , we prepared another code generator GENCODELBN similar to GENCODEN by changing projection operators.The 47 integrals for ∆LB 8 were then numerically evaluated.The numerical values of all finite integrals in (99) are listed in Table 3.Many of them are borrowed from [40].The tenth-order finite magnetic moment amplitude is updated in this work as in units of (α/π) 5 .The factor two of each time-reversal-symmetric diagram is included in the numerical value of ∆M G .The improvement of the numerical value of ∆M 10 leads to a new tenth-order contribution A (10) 1 of (8)., and ∆dm 2n are the sum of the finite magnetic moment amplitudes, the sum of the finite parts of vertex and wave-function renormalization constants, and the sum of the finite parts of the mass-renormalization constants, respectively, all derived from the 2n th -order diagrams without a fermion loop of the QED perturbation theory.They are given in units of (α/π) 2n .∆M 10 is newly calculated for this paper.∆M 8 is derived from the near-exact result Equation ( 5) of [26].

Conclusions
An overview of the current status of the standard model prediction of the electron anomalous magnetic moment a e is given.The fine-structure constant determined from measurements in atomic physics is the dominant source of the uncertainty of a e .Both numerical and analytic works on the QED contribution have succeeded in reducing its uncertainty.By now, the hadronic contribution is the second largest source of uncertainty in the standard model prediction of a e .The method of numerical computation of the higher order QED contributions to a e , especially those from diagrams without a fermion loop, is described in some detail.

Figure 1 .
Figure 1.Second-order vertex diagram.There is only one diagram.The straight and wavy lines represent electron and photon propagators, respectively.Reprinted from [12].

Figure 2 .
Figure2.Fourth-order vertex diagrams.There are seven diagrams in total.The time-reversal diagrams of (a,c) are not shown.The solid and wavy lines represent electron and photon propagators, respectively.Reprinted from[12].

Figure 3 .
Figure 3. Sixth-order vertex diagrams.There are 72 diagrams in total, and they are divided into five gauge-invariant sets.Typical diagrams from each set are shown as (a-e).There are (a) 3 diagrams, (b) 1 diagram, (c) 12 diagrams, (d) 6 diagrams, and (e) 50 diagrams.The solid and wavy lines represent electron and photon propagators, respectively.Reprinted from [12].

Figure 6 .
Figure 6.Fourth-order self-energy-like diagrams, M 4a and M 4b .The solid and wavy lines represent electron and photon propagators, respectively.

Figure 7 .
Figure7.A tenth-order self-energy-like diagram.This is one of the diagrams V of Figure5.The solid and wavy lines represent electron and photon propagators, respectively.Reprinted from[39].
2. inclusion:S α ⊂ S β (a)If S α is a magnetic moment, S β cannot be a self-mass: I S α I S β .(b)Both S's are self-masses: R S β /S α R α .

Table 2 .
Standard model contributions to the electron anomalous magnetic moment a e in units of 10 −12 .