Radiative-nonrecoil corrections of order alpha^2 (Z alpha) E_F to the hyperfine splitting of muonium

We present results for the corrections of order alpha^2 (Z alpha) E_F to the hyperfine splitting of muonium. We compute all the contributing Feynman diagrams in dimensional regularization and a general covariant gauge using a mixture of analytical and numerical methods. We improve the precision of previous results.


I. INTRODUCTION
Muonium is the hydrogenlike bound state of a positive muon and an electron. Unlike hydrogen, or any other bound state involving hadrons, muonium is free from the complications introduced by the finite size or the internal structure of any of its constituents. Therefore, it allows for a very precise test of bound-state QED, and can be used to restrict models of physics beyond the Standard Model. Measurements of the ground-state hyperfine splitting of muonium are used to extract the muon to electron mass ratio m µ /m e and the muon to proton magnetic moment ratio µ µ /µ p [1]. The value of µ µ /µ p is required for obtaining the muon anomalous magnetic moment from experiment [2]. In addition, the hyperfine splitting can also be used to determine the fine structure constant α. For a review of the present status and recent developments in the theory of light hydrogenic atoms, see [3,4].
The leading-order hyperfine splitting is given by the Fermi energy E F (defined in Eq. (1)). Its corrections are organized as a perturbative expansion in powers of three parameters: Zα, describing effects due to the binding of an electron to a nucleus of atomic number Z; α (frequently accompanied by 1/π) from electron and photon self-interactions; and the ratio of electron to nucleus masses, m/M . The main theoretical uncertainty comes from three types of yet unknown corrections: single-logarithmic and non-logarithmic corrections of order α(Zα) 2 (m/M )E F , and non-logarithmic corrections of order α 2 (Zα)(m/M )E F and (Zα) 3 (m/M )E F (some terms are known for the first case [5]).
In this paper we focus on the second-order radiativenonrecoil corrections to the hyperfine splitting (of order α 2 (Zα)E F ). The total result for these corrections was found by Eides and Shelyuto [6] and Kinoshita and Nio [7]. Our result improves their precision by over an order of magnitude. Our central value is slightly lower than, but compatible with, that of [6].
In Sec. II we present the details of our approach, and in Sec. III we present our results. In Appendix A we show analytic results for two sets of diagrams.

II. EVALUATION
We consider an electron of mass m orbiting a nucleus of mass M and atomic number Z. In this paper we consider the nucleus to be a muon, but we will keep Z explicit in order to distinguish between the binding contributions (Zα) and the radiative ones (α).
We are interested in corrections to the hyperfine splitting of the ground state of muonium of order α 2 (Zα)E F and leading order in m/M , where Here g is the gyromagnetic factor of the nucleus 1 (in our case, a muon, but our final result in Eq. (19) applies to any hydrogenlike atom). In order to compute these corrections, we consider the scattering amplitude where u(p) is the spinor for the electron, v(P ) is the spinor for the muon, p = (m, 0) and P = (M, 0). Q 1 and Q 2 are given by the Feynman rules describing the sum of the diagrams shown in Figs. 1 and 2. In these figures, the sum of the direct and crossed interactions between the electron and the muon is represented by a dotted line, as shown in Fig. 3. We define a bound-state wave function ψ = uv, so that Eq. (2) becomes I II III IV V Figure 1: The different sets of vacuum polarization diagrams and light-by-light diagrams (set IV). Each set represents the drawn diagram plus all the possible permutations of its pieces.
Depending on the relative alignment of the spins of the constituent particles, an S state can either belong to the J = 1 triplet or the J = 0 singlet. The triplet and singlet states are often denoted by the prefixes ortho-and para-, respectively, and their wave functions are given by [8] where ξ is the polarization vector. We average over the directions of ξ by considering the four-vector ξ ≡ (0, ξ) and using the identity We use dimensional regularization with d = 4−2ǫ dimensions. Thus, an important issue is the definition of γ 5 , which is an intrinsically four-dimensional object. Since we do not have to evaluate traces with an odd number of γ 5 matrices, we can treat them as anticommuting. The energy shift created by the radiative corrections depicted in Figs. 1 and 2, for either the singlet or triplet configurations, is given by where |ψ n (0)| 2 = (Zαµ) 3 /(πn 3 ) is the squared modulus of the wave function of a bound S state with principal quantum number n and reduced mass µ. The hyperfine splitting is then simply In order to evaluate the loop integrals represented by the Feynman diagrams we use the method of regions [9] to construct an expansion in the small ratio m/M . There are several possible contributing regions, where one or more of the loop momenta scale like m or M . However, we are only interested in the leading order in m/M , which is given by the region where all loop momenta scale like m. If k ∼ m, we can expand the contribution from the muon line in the sum of the direct and crossed diagrams of Fig. 3,  where We used the equation of motion to set some terms in the numerator to zero, and we arranged the terms in the expansion in such a way that the three different Dirac structures that are important for the calculation of the hyperfine splitting appear explicitly. We will now see that only T 2 can contribute to the splitting. Consider the Dirac structure of ψT 1 and anticommute the gamma matrices, for both para and ortho states: Now we can write Using the expressions in Eqs. (13) and (14) it is easy to see that χ para T1 ψ † para = χ ortho T1 ψ † ortho (after averaging over polarizations). This means that T 1 gives the same contribution for para and ortho states. Therefore, when we subtract these contributions in order to compute the hyperfine splitting, they cancel out.
If we consider T 2 instead, defining χ T2 in analogy with Eqs. (13) and (14), we can see that χ para T2 ψ † para = χ ortho T2 ψ † ortho , so this term will not cancel in the subtraction. The difference between the para and ortho states comes solely from terms in χ ortho T2 that are totally antisymmetric in α and β. Therefore, when we consider the Dirac structure of T 3 , which is but a symmetrization of that of T 2 , these terms will vanish, and so T 3 will give no contribution to the hyperfine splitting either.
Thus, we have seen that the only term that contributes to the hyperfine splitting is (16) This is valid at all orders of alpha, and all orders in m/M . We can then substitute the scalar part of the nucleon propagator by a Dirac delta in all our calculations, since we are only interested in the leading order in m/M and We used dimensional regularization, and renormalized our results using the on-shell renormalization scheme. For all the photon propagators in Figs. 1 and 2 we used a general covariant R ξ gauge. The overall cancellation of the dependence on the gauge parameter in the final result provides us with a good check for our calculations.
We used the program qgraf [10] to generate all of the diagrams, and the packages q2e and exp [11,12] to express them as a series of vertices and propagators that can be read by the FORM [13] package MATAD 3 [14]. Finally, MATAD 3 was used to represent the diagrams in terms of a set of scalar integrals using custom-made routines. In this way, we represented the amplitude M in terms of several thousand different scalar integrals. These integrals can be expressed in terms of a few master integrals by means of integration-by-parts (IBP) identities [15]. We used the so-called Laporta algorithm [16,17] as implemented in the Mathematica package FIRE [18], to reduce the problem to 32 master integrals. The master integrals for this calculation are the same ones we found in [19]. All definitions and results for the integrals can be found in this  Fig. 1 and those of [22][23][24]. Numbers ending in an ellipsis indicate an analytic result, which we show in Appendix A. No error was given for the numerical result of set I in [23]. reference. However, one change was made for this calculation. In order to obtain better numerical precision, we performed a change of basis, so that instead of working with I 14 = F (1, 0, 0, 0, 1, 1, 1, 1) we worked with which was obtained using the Mathematica package FIESTA 1.2.1 [20] with integrators from the CUBA library [21].

Our final result for the hyperfine splitting is
This correction was also found by Eides and Shelyuto [6], and by Kinoshita and Nio [7]. Their results are δE hf s = 0.7679(79) · α 2 (Zα) Our result is a little over one order of magnitude more precise than that of [6], and almost three orders of magnitude more precise than the one in [7]. Our central value is slightly lower than in [6], by about 1.5σ (taking as σ the larger error). It agrees with the result of [7] within its much larger error estimate. For the ground state of muonium, our result reads We compared our results for the individual diagrams and those found in the literature [6,7,[22][23][24]. Our results for the gauge-invariant sets of diagrams of Fig. 1 are presented in Table I. For the diagrams of Fig. 2 we chose the Fried-Yennie gauge [25,26], in which all diagrams are infrared finite. Our results are presented in Table II.
The sum of all central values in the second column of Tables I and II gives the coefficient 0.77099 in Eq. (19). The error of that result is however not obtained from the sum of the errors of the diagrams in the tables. Once we decompose the problem into the calculation of master integrals, the diagrams are no longer independent, as the same master integral contributes to several different diagrams. Thus, to find the error of our total result, we first sum all diagrams and then sum all the errors of the integrals in quadrature. We found new analytic results for diagrams g and f , shown in Table II. For completeness, the known analytic results for sets II and III of the vacuum polarization diagrams are given in Appendix A as well.
We found no discrepancies between our results for the diagrams of Fig. 1 and the ones of [22][23][24], but we found significant differences in the rest of the diagrams between the results of [6] and ours. They affect all diagrams except diagrams a, e, h, l, o, and q. The biggest discrepancies are in diagrams b and c, and they are of the order of 30σ. However, most of the differences cancel when summing the diagrams. In particular, there are almost exact cancellations between the differences in diagrams b and c, k and l, and n and p.
The reason for the discrepancies (and their cancellations) is most likely the different treatment of infrared divergences in [6] and this paper. In [6], the Fried-Yennie gauge was set from the beginning, and all spurious infrared divergences were canceled before the integration over the diagram's loop momenta, which was performed in four dimensions. In our calculation, we used a general gauge parameter, and dimensional regularization to deal with infrared divergences, which would only vanish after setting the gauge in the final expression. As noted in [27], there is a difference between setting the Fried-Yennie gauge and sending the infrared regulator to zero before or after integration. It is not surprising then that we obtained different results than [6] for gauge-dependent diagrams, but that most of the differences cancel in the final, gauge-invariant result, making it compatible with the previous calculation.
Using the setup of the calculation of the hyperfine splitting one can also find the Lamb shift, as it is given by We obtained in this way the same results as in [19] 2 .