The $O(\alpha^2)$ Initial State QED Corrections to $e^+e^-$ Annihilation to a Neutral Vector Boson Revisited

We calculate the non-singlet, the pure singlet contribution, and their interference term, at $O(\alpha^2)$ due to electron-pair initial state radiation to $e^+ e^-$ annihilation into a neutral vector boson in a direct analytic computation without any approximation. The correction is represented in terms of iterated incomplete elliptic integrals. Performing the limit $s \gg m_e^2$ we find discrepancies with the earlier results of Ref.~\cite{Berends:1987ab} and confirm results obtained in Ref.~\cite{Blumlein:2011mi} where the effective method of massive operator matrix elements has been used, which works for all but the power corrections in $m^2/s$. In this way, we also confirm the validity of the factorization of massive partons in the Drell-Yan process. We also add non-logarithmic terms at $O(\alpha^2)$ which have not been considered in \cite{Berends:1987ab}. The corrections are of central importance for precision analyzes in $e^+e^-$ annihilation into $\gamma^*/Z^*$ at high luminosity.

The initial state QED corrections (ISR) to e + e − annihilation are of crucial importance for the experimental analyzes at LEP [3] and for planned projects like the ILC and CLIC [4], the FCC ee [5], muon colliders [6], and notably also at e + e − Higgs factories using the process e + e − → ZH 0 , as well as e + e − → tt. The initial state corrections have been carried out analytically to O (αL) 5 in the leading logarithmic series using the structure function method [7]. A small z-resummation has been performed in [8]. In Ref. [1] the O(α 2 ) corrections were calculated neglecting terms of O m 2 s ln( m 2 s ) . Here s is the cms energy squared and m = m e . These corrections are used in analysis codes such as TOPAZ0 [9] and ZFITTER [10]. The initial state QED corrections can be written in terms of the following functions which yield the respective differential cross sections by with σ e + e − (s ) the scattering cross section without the ISR corrections, α ≡ α(s) the fine structure constant and z = s /s, where s is the invariant mass of the produced (off-shell) γ/Z boson. In Ref. [2] the method of massive operator matrix elements (OMEs) [11] has been applied to calculate the O(α 2 ) corrections, factorizing the process into universal massive contributions and the massless Wilson coefficients of the Drell-Yan process [12,13]. This method has been known to work in the case of external massless fields, including the non-logarithmic contributions in QCD, cf. [11,14]. However, the results of Refs. [1] and [2] were found to disagree.
One way to find out the correct answer is to perform the direct analytic calculation of the corresponding contributions without doing any approximation. This is done in the present paper for three of the subprocesses of the O(α 2 ) corrections related to fermion-pair production, i.e. for the processes II-IV of Ref. [1]. Like in Refs. [1,2] we will consider the case of pure vector couplings in the following. It is known already from the massless case [12] that in some of the processes axialvector-vector terms receive different corrections if compared to the pure vector or axialvector case. These aspects will be presented in Ref. [15]. In the calculation we used the packages FORM, Sigma, HarmonicSums, HolonomicFunctions [17][18][19][20][21][22][23][24][25][26][27][28][29] and private implementations [30]. The complete results have an iterative integral representation. We compare the exact result numerically with the one obtained in the limit ρ = m 2 /s → 0. Both results agree better than a relative deviation of 10 −7 at s = M 2 Z , as expected by neglecting the power corrections. The result is given in terms of the variable z and powers of the logarithm L = ln(s/m 2 ). The logarithmic corrections in [1,2] agree.
In the calculation, the phase space integrals can be mapped to fourfold scalar integrals. Here one depends on s /s, with s the mass squared of the emitted fermion pair. Unlike the case of s and s , the latter invariant is not large against m 2 everywhere. Three of the four integrals can be consecutively obtained yielding results which contain different functions whose arguments involve square-roots. It is then useful to construct a basis of the contributing root-valued letters and to perform the last integral over it, using differential field methods [27].
Here also nested roots have to be transformed to single roots before. By this the complete integrals are at most triple iterated over an alphabet including also new types of square-root letters. Finally, the individual terms are regularized such that the expressions can be expanded in the ratio ρ term by term [16]. Typical letters are Here y denotes the next integration variable. After having performed suitable regularizations, one may also expand in ρ before performing the last integral. Processes II and III have been also considered primarily for massless quark-antiquark pair production in Refs. [31,32] in 4-dimensional calculations. Here a quark mass serves as a regulator since the massless limit is aimed at from the beginning. Neglecting the mass terms not needed to regularize the corresponding integrals leads to simpler integrands, which finally integrate to polylogarithms directly and the logarithmic contributions in L are obtained correctly. In the case discussed at present, however, the finite electron mass is physical and the expansion in m 2 is only possible if all terms which contribute to the final result are retained. This may require a deeper expansion in m 2 than the one performed in Ref. [1]. Terms which can be safely neglected are of O ρ ln k (ρ) , k = 0, 1, 2 in the result.
For process II in [1] we find the difference term Here the original lower bound of the integral 4m 2 e /(s(1 − √ z) 2 ) can be set to zero since its contribution is of the order of the neglected terms. From this integral it follows that terms containing higher powers in 1/(1 − z) are appearing, which were not contained in [1,31]. They emerge from further terms in the mass expansion that cannot be neglected.
In the case of process III, Ref. [1], takes the results from [32] in which, however, mass terms being necessary here, were neglected beforehand since the result concerned a massless quark calculation [33]. Note also that the pure singlet interference term [32] was taken with the wrong sign.
The difference terms, δ i , i = II, III, IV, for e + e − pair emission between the present results and those of Ref. [1], after the analytic expansion of the complete expression in m 2 /s including the constant term, read : with δ PS interf the term given in (B.22) of [12] setting C F = T F = 1, α s = α. Here, Li k (z) and S p,k (z) denote the polylogarithm and Nielsen integrals, respectively, cf. [35]. We remark that we agree with the result on the interference terms in the pure singlet case [32] which has been found there already to be regularization scheme invariant; it also agrees with the massless result [12]. Numerically, it turns out that the deviation due to the δ PS interf is small against the other differences for the pure singlet term. We agree with the result for µ + µ − emission for the non-singlet term (process II) of Refs. [1,31], see also Ref. [15]. This contribution has also a representation using massive OMEs. Here, however, the external fermion lines are massless since m e m µ , see Ref. [36]. The terms of O(1/(1 + z) 2(3) ) given in [1] for process IV have not been observed in the respective contribution of the OMEs in [2]. Terms of this kind are only expected for (anti)particle-(anti)particle scattering.
The relative deviations for the results for processes II-IV in the present calculation and [1] are shown in Figure 1, where ∆ (2) denotes the ratio of δ i and the corresponding complete O(α 2 ) correction for i = II, III, IV. All illustrations are made for z < 1. The relative differences reach from +25 to -60% for z ∈ [10 −5 , 1]. Here we have changed the term ln(z)/(1−z) 2 → ln 2 (z)/(1−z) 2 in Eq. (2.43) in [1] which appears twice (suggesting a typo), such that this term is only logarithmic but not linear divergent for z → 1 and thus integrable. Otherwise the difference would be even larger. For the non-logarithmic terms, not all contributions have been considered in [1]. These are the graphs B and their interference terms with the non-singlet (A) and pure singlet terms (C and D) in [12]. We have recalculated them in the massive case. Note that the interference term between the graphs A and B only contributes to the axialvector term. For these processes there is no massive OME. Since the massive OMEs contain all massive corrections in the limit m 2 /s → 0, the only contributions are from the massless Wilson coefficients. We find in the corresponding massive calculation the massless results given in [12], which is a further confirmation of the formalism presented in Ref. [11] including the constant terms. This has also been observed e.g. in the case of the massive asymptotic two-loop corrections to the deep-inelastic structure function F L (x, Q 2 ) [11,37]. Figure 2 shows the different contributions at O(α 2 ) of initial state e + e − pair production to γ * /Z * -boson production. The dominant contributions come from the pure singlet and non-singlet terms; other contributions are smaller but not negligible at the 0.1% level in the radiator function. For large values of z = s /s the non-singlet terms are dominant, whereas for z ∼ < 0.03 the pure singlet contributions dominate. Analogous contributions to those considered here have been calculated at three-loop order for the massive OMEs in QCD with external massless parton lines in [36]. When Ref. [2] was published, we could not explain the differences to Ref. [1] and we tended to assume that the OME method might have a problem in case of massive external states, which is now proven not to be the case. The factorization of massive initial states for the Drell-Yan process [38] is also observed in the case discussed here.

ACKNOWLEDGMENTS
This paper is dedicated to the memory of our colleague W.L. van Neerven. We would like to thank J.C. Collins, J.H. Kühn, G. Passarino, C. Schneider and G. Sterman for discussions. This project has received funding from the European Unions Horizon 2020 research and innovation programme under the Marie Sklodowska-Curie grant agreement No. 764850, SAGEX, and COST action CA16201: Unraveling new physics at the LHC through the precision frontier and from the Austrian FWF grants P 27229 and P 31952 in part.