A Evaluating the G-functions

## Abstract

We present the two-mass QCD contributions to the pure singlet operator matrix element at three loop order in -space. These terms are relevant for calculating the structure function at as well as for the matching relations in the variable flavor number scheme and the heavy quark distribution functions at the same order. The result for the operator matrix element is given in terms of generalized iterated integrals that include square root letters in the alphabet, depending also on the mass ratio through the main argument. Numerical results are presented.

DESY 17–195

DO–TH 17/32

November 2017

The two-mass contribution to the three-loop

pure singlet operator matrix element

J. Ablinger, J. Blümlein, A. De Freitas, C. Schneider,

and K. Schönwald

Research Institute for Symbolic Computation (RISC),

Johannes Kepler University, Altenbergerstraße 69, A–4040, Linz, Austria

Deutsches Elektronen–Synchrotron, DESY,

Platanenallee 6, D-15738 Zeuthen, Germany

## 1 Introduction

Starting at 2-loop order, massive operator matrix elements (OME) , which are the transition matrix elements in the variable flavor number scheme (VFNS), receive two-mass contributions [1, 2]. This also applies to the Wilson coefficients in deeply inelastic scattering. The single mass contributions to the pure singlet (PS) OME, , have been calculated in Ref. [3]. In the present paper we present the corresponding 2-mass contributions. They occur first at 3–loop order. Previously, we have calculated already the fixed moments of this OME for the Mellin variable , to , with the mass ratio of the heavy quarks squared in [4, 2] using the package Q2E [5, 6], as well as the two-mass contributions to the 3-loop OMEs in the flavor non-singlet cases and and for , cf. [2], both in - and in -space. The 2-mass contributions to are in preparation [7]. Various other 3-loop single mass contributions have also been completed, cf. [13, 14, 16, 15, 8, 9, 10, 12, 11, 17] and all the logarithmic contributions are known to this order [18]. Furthermore, for the OME , all diagrams consisting of contributions that can be obtained in terms of first order factorizing differential or difference equations have been calculated [20, 19, 21]. For all OMEs a series of moments has been calculated in the single mass case in Ref. [22] using the code MATAD [23]. From the single pole terms of all the OMEs, we have derived the contributing 3-loop anomalous dimensions, cf. e.g. [3, 24].

We perform the calculation of the 2-mass part of the OME mainly in -space1, using only some elements of the -space formalism. In the presnt case one obtains first order factorizable expressions in -space, but not in -space. This implies that the -space solution cannot be given by sum and product expressions only.

The paper is organized as follows. In Section 2 we present the renormalized pure-singlet OME in the 2-mass case. Details of the calculation are given in Section 3. In Section 4 the result of the calculation is presented and numerical results are given in Section 5. Section 6 contains the conclusions. Details on new iterated integrals emerging in the representation, the -functions, and a series of fixed moments as a function of the mass ratio of the two heavy quarks for the unrenormalized OME are given in the Appendix.

## 2 The renormalized 2-mass pure singlet OME

The generic pole structure for the PS three–loop two–mass contribution is given by [2]

 ^^~A(3),PSQq = 83ε3γ(0)gq^γ(0)qgβ0,Q+1ε2[2γ(0)gq^γ(0)qgβ0,Q(L1+L2)+16^γ(0)qg^γ(1)gq−43β0,Q^γPS,(1)qq] (2.1) +1ε[γ(0)gq^γ(0)qgβ0,Q(L21+L1L2+L22)+{18^γ(0)qg^γ(1)gq−β0,Q^γPS,(1)qq}(L2+L1) +13^~γ(2),PSqq−8a(2),PSQqβ0,Q+^γ(0)qga(2)gq]+~a(3),PSQq(m21,m22,μ2),

where we used the short hand notation2

 ^γij = γij(NF+2)−γij(NF) (2.2) ^~γij = γij(NF+2)NF+2−γij(NF)NF. (2.3)

The tilde in indicates that we are considering only the genuine two-mass contributions, and the double hat is used to denote a completely unrenormalized OME. Here the ’s are anomalous dimensions at loops, , and

 L1=ln(m21μ2),L2=ln(m22μ2), (2.4)

where and are the masses of the heavy quarks, and is the renormalization scale. Our goal is to compute the term .

In the –scheme, renormalizing the heavy masses on-shell, the renormalized expression is given by

 ~A(3),¯¯¯¯¯¯¯MS,PSQq = −γ(0)gq^γ(0)qgβ0,Q(14L22L1+14L21L2+13L31+13L32) (2.5) +{−116^γ(0)qg^γ(1)gq+12β0,Q^γPS,(1)qq}(L22+L21) +{4a(2),PSQqβ0,Q−12^γ(0)qga(2)gq−14β0,Qζ2γ(0)gq^γ(0)qg}(L1+L2) +8¯¯¯a(2),PSQqβ0,Q−^γ(0)qg¯¯¯a(2)gq+~a(3),PSQq(m21,m22,μ2).

The transition relations for the renormalization of the heavy quarks in the -scheme is given in [3], Eq. (5.100), but it only applies to the equal mass case since for the unequal mass case the first contributions emerge at 3–loop order. In Eqs. (2.1) and (2.5), and represent the terms of the two-loop OMEs and , respectively, while and represent the corresponding terms, cf. Refs. [26, 27, 28, 29, 30]. Here and in what follows, denotes the Riemann -function.

## 3 Details of the calculation

### 3.1 The basic formalism

There are sixteen irreducible diagrams contributing to , which are shown in Figure 1. The unrenormalized operator matrix element is obtained by adding all the diagrams and applying the quarkonic projector to the corresponding Green function ,

 ^^~A(3)Qq(m21μ2,m22μ2,ε,N)=Pq^GijQ≡δij4Nc(Δ.p)−NTr[⧸p^GijQ], (3.1)

where is the momentum of the on-shell external massless quark (), is a light-like -vector, with , and are the color indices of each external leg, and is the number of colors. The diagrams, , are calculated directly within dimensional regularization in dimensions. These diagrams have the following structure,

 Di(m1,m2,N) = (m2μ2)32ε~Di(η,ε,N) (3.2) = (m2μ2)32ε∫10dxxN−1^Di(η,ε,x), (3.3)

where is the Mellin variable appearing in the Feynman rules for the operator insertions, cf. [22, 2], is the mass of the heavy quark where the operator insertion is not sitting, and

 η=m22m21, (3.4)

with , i.e. .

In previous publications where single-mass OMEs have been computed [8, 18, 9, 10, 3, 15, 16, 20, 12], and in our recently published two-mass calculations [2], we have always given the results both in - and -space. In the case of the two-mass pure singlet OME, finding a general -space result at three loops, turns out to be rather cumbersome. We will, therefore, present our result only in -space, see Eq. (3.3), which is anyway all we need in order to obtain the corresponding contribution to the structure function for large values of , as well as the contribution to the variable flavor number scheme. In most of the applications one finally works in -space.

The diagrams on the first row of Figure 1 all give the same result, and the diagrams on the second row are related to the diagrams on the first row by the exchange , i.e.,

 Di(m1,m2,N) = D1(m1,m2,N)fori=2,3,4. (3.5) Di(m1,m2,N) = D1(m2,m1,N)fori=5,6,7,8. (3.6)

All of these diagrams vanish for odd values of .

A relation similar to (3.6) holds between the diagrams on the third row of Figure 1 and those on the fourth row. Furthermore, the diagrams on these two last rows which differ only in the direction of a fermion arrow are related by a factor of relative to each other. Specifically one has

 D10(m1,m2,N) = D9(m1,m2,N). (3.7) Di(m1,m2,N) = (−1)ND9(m1,m2,N)fori=11,12. (3.8) Di(m1,m2,N) = D9(m2,m1,N)fori=13,14. (3.9) Di(m1,m2,N) = (−1)ND9(m2,m1,N)fori=15,16. (3.10)

We can therefore write the whole unrenormalized pure singlet operator matrix element solely in terms of diagrams 1 and 9:

 A(3),PSQq(N) = 4D1(m1,m2,N)+2(1+(−1)N)D9(m1,m2,N) (3.11) +4D1(m2,m1,N)+2(1+(−1)N)D9(m2,m1,N).

In Eq. (3.11) the use of the projection (3.1) is implicit.

All of the diagrams contain a massive fermion loop with an operator insertion (Figures 2 and 2) and a massive bubble without the operator (Figure 2).

The latter can be rendered effectively massless by using a Mellin–Barnes integral [31, 32, 33, 34, 35]

 Iμν,aba1(k) = g2sTF4π(4π)−ε/2(kμkν−k2gμν) (3.12) ×∫+i∞−i∞dσ(m2μ2)σ(−k2)ε/2−σΓ(σ−ε/2)Γ2(2−σ+ε/2)Γ(−σ)Γ(4−2σ+ε),

where and are the respective Lorentz indices of the external legs, and the color indices, is the external momentum, is the mass of the fermion, which can be either or , is the strong coupling constant, and in , with the number of colors.

Diagram 1 can then be calculated using the following expression for the massive fermion bubble containing the vertex operator insertion (Figure 2) :

 Iμν,abb1(k) = 16δabTFg2s(Δ.k)N−2(4π)D/2Γ(2−D/2)∫10dxxN(1−x)(Δ.k)Δμkν−k2ΔμΔν(m2−x(1−x)k2)2−D/2.

Likewise, diagram 9 can be calculated using the following expression for the bubble containing the operator insertion on the fermion line (Figure 2):

 Iμν,abb2(k) = 4δabTFg2s(Δ.k)N−2(4π)D/2∫10dxxN−2(1−x)[ (3.14) −2(x(1−x)(gμνk2−2kμkν)+m2gμν)x2Γ(3−D/2)(Δ.k)2(m2−x(1−x)k2)3−D/2 +Γ(2−D/2)(2Nx+1−N)x(kμΔν+kνΔμ)(Δ.k)(m2−x(1−x)k2)2−D/2 +Γ(2−D/2)((N−1)(1−2x)−Dx)xgμν(Δ.k)2(m2−x(1−x)k2)2−D/2 −Γ(1−D/2)N−11−x(N(1−x)−1)ΔμΔν(m2−x(1−x)k2)1−D/2].

After inserting these expressions, we introduced the corresponding projector for a quarkonic OME given in Eq. (3.1), and performed the Dirac matrix algebra and the trace arising in the numerator of the diagrams using the program FORM [36]. This leads to a linear combination of integrals, the denominators of which were then combined using Feynman parameters. The momentum integrals were then performed, and we obtained expressions where one of the Feynman parameter integrals is already in the form of a Mellin transform. Therefore, we left this Feynman parameter un-integrated, and integrated the remaining ones, after which we obtained the following expression for diagram 1,

 D1(m1,m2,N)=−128CFT2F(1+(−1)N)(J1−J2), (3.15)

with and , and

 J1 = (m21μ2)32εΓ(N−1)Γ(N+ε2)∫10dxxN+ε2(1−x)1+ε2B1(ηx(1−x)), (3.16) J2 = (m21μ2)32εΓ(N)Γ(N+1+ε2)∫10dxxN+ε2(1−x)1+ε2B1(ηx(1−x)). (3.17)

Here is the following contour integral,

 B1(ξ)=12πi∫i∞−i∞dσξσΓ(−σ)Γ(−σ+ε)Γ(σ−3ε2)Γ(σ−ε2)Γ2(σ+2−ε)Γ(2σ+4−2ε). (3.18)

For diagram 9, we get

 D9(m1,m2,N) = 64CFT2F{−14(ε+2)[−2J3+2ηJ4+(2N+2+ε)J5 (3.19) −(N−1)J6]+N(N−1)(J7−J8)−(N−1)(J9−J10)},

where

 J3 = (m21μ2)32εΓ(N+1)Γ(N+2+ε2)∫10dxxN+ε2(1−x)1+ε2B2(ηx(1−x)), (3.20) J4 = (m21μ2)32εΓ(N+1)Γ(N+2+ε2)∫10dxxN−1+ε2(1−x)ε2B3(ηx(1−x)), (3.21) J5 = (m21μ2)32εΓ(N+1)Γ(N+2+ε2)∫10dxxN+ε2(1−x)1+ε2B1(ηx(1−x)), (3.22) J6 = (m21μ2)32εΓ(N+1)Γ(N+2+ε2)∫10dxxN−1+ε2(1−x)1+ε2B1(ηx(1−x)), (3.23) J7 = (m21μ2)32εΓ(N−1)Γ(N+ε2)∫10dxxN−1+ε2(1−x)2+ε2B4(ηx(1−x)), (3.24) J8 = (m21μ2)32εΓ(N)Γ(N+1+ε2)∫10dxxN−1+ε2(1−x)2+ε2B4(ηx(1−x)), (3.25) J9 = (m21μ2)32εΓ(N−1)Γ(N+ε2)∫10dxxN−1+ε2(1−x)1+ε2B4(ηx(1−x)), (3.26) J10 = (m21μ2)32εΓ(N)Γ(N+1+ε2)∫10dxxN−1+ε2(1−x)1+ε2B4(ηx(1−x)), (3.27)

with

 B2(ξ) = 12πi∫i∞−i∞dσξσΓ(−σ)Γ(−σ+ε)Γ(σ−3ε2)Γ(σ+1−ε2)Γ2(σ+2−ε)Γ(2σ+4−2ε), (3.28) B3(ξ) = 12πi∫i∞−i∞dσξσΓ(−σ)Γ(−σ−1+ε)Γ(σ+1−3ε2)Γ(σ+1−ε2) (3.29) 12πi∫i∞−i∞dσ×Γ2(σ+3−ε)Γ(2σ+6−2ε), B4(ξ) = 12πi∫i∞−i∞dσξσΓ(−σ)Γ(−σ+ε)Γ(σ−3ε2)Γ(σ−1−ε2)Γ2(σ+2−ε)Γ(2σ+4−2ε).

We want to compute the diagrams, and therefore the integrals (3.18, 3.28LABEL:B4), as an expansion in up to . Notice that there is always a factor consisting on a ratio of –functions depending on in these integrals, and there are also additional factors of in Eq. (3.19). After the aforementioned expansion is performed, some of the integrals will be left with a factor of the form

 1N+l,withl∈{−1,0,1}. (3.31)

In order to get our results really in terms of a Mellin transform, these factors need to be absorbed into the integral in , which can be accomplished using integration by parts,

We will see later how this is manifested in the final result.

### 3.2 Computation of the contour integrals

We have seen that the diagrams where the operator insertion lies on the heaviest quark are related to the ones where the insertion lies on the lightest internal quark by the exchange , which means that in order to go from the former diagrams to the latter, we also need to do the change . Therefore, while in the case of diagrams and , we need to evaluate the contour integrals (3.18, 3.28LABEL:B4) with

 ξ=ηx(1−x), (3.34)

in the case of diagrams and , we have instead

 ξ=1ηx(1−x). (3.35)

As we will see, the contour integrals arising in the case of Eq. (3.34) will require a different treatment to those in the case of Eq. (3.35).

The integrals (3.18, 3.28LABEL:B4) can be computed with the help of the Mathematica package MB [37], together with the add-on package MBresolve [38]. These packages allow us to resolve the singularity structure in of these integrals by taking residues in , after which we can perform the expansion in . This leads to expressions of the form

 Bi(ξ)=B(ε)i(ξ)+B(0)i(ξ),i=1,2,3,4, (3.36)

where the functions are the sum of residues. For example, in the case of , we must take residues at and , which leads to

 B(ε)1(ξ) = ξε/2Γ(−ε)Γ(−ε2)Γ(ε2)Γ2(2−ε2)Γ(4−ε) (3.37) +ξ3ε/2Γ(−3ε2)Γ(−ε2)Γ(ε)Γ2(2+ε2)Γ(4+ε).

Once the residues are taken, it is possible to expand the original integrals in . This leads to the contour integrals that we denote by , for which the MBresolve package finds a contour consisting of a line parallel to the imaginary axis with a real value for 3. The functions are given by

 B(0)1(ξ) = 12πi∫−1/2+i∞−1/2−i∞dσξσΓ2(−σ)Γ2(σ)Γ2(σ+2)Γ(2σ+4)+O(ε), (3.38) B(0)2(ξ) = 12πi∫−1/2+i∞−1/2−i∞dσξσΓ2(−σ)Γ(σ)Γ(σ+1)Γ2(σ+2)Γ(2σ+4)+O(ε), (3.39) B(0)3(ξ) = 12πi∫−1/2+i∞−1/2−i∞dσξσΓ(−σ)Γ(−σ−1)Γ2(σ+1)Γ2(σ+3)Γ(2σ+6)+O(ε), (3.40) B(0)4(ξ) = 12πi∫−1/2+i∞−1/2−i∞dσξσΓ2(−σ)Γ(σ)Γ(σ−1)Γ2(σ+2)Γ(2σ+4)+O(ε). (3.41)

In order to calculate the contour integrals (3.383.41), we need to close the contour either to the right or to the left, depending on the convergence of the resulting sum of residues. The integrand always contains a -function in the denominator that can be rewritten using Legendre’s duplication formula,

 Missing or unrecognized delimiter for \left (3.42)

which means that if we close the contour to the right, the sum of residues will be convergent if and only if , while if we close the contour to the left, the sum of residues will be convergent if and only if . Therefore, in the case where , we just need to close the contour to the left, since , and therefore . The calculation in the case where , on the other hand, is a bit more complicated, since in this case we need to split the integration range in into the different regions where can be bigger or smaller than 4,

 ηx(1−x) < 4,forx∈(η−,η+), (3.43) ηx(1−x) > 4,forx∈(0,η−)orx∈(η+,1), (3.44)

where

 η±=12(1±√1−η). (3.45)

Therefore, in this case, we will need to close the contour to the right when , and to the left when or .

Closing the contour to the right and summing residues, we obtain

 B(0)1(ξ) = −(ζ26+1427)ln(ξ)−136ln3(ξ)+536ln2(ξ)+518ζ2−ζ33+8281 (3.46) +∞∑k=1ξkΓ2(k+2)k2Γ(2k+4)(2S1(2k+3)−2S1(k)+2k(k+1)−ln(ξ)), B(0)2(ξ) = −112ln2(ξ)+518ln(ξ)−ζ26−1427+∞∑k=1ξkΓ2(k+2)kΓ(2k+4)(2S1(2k+3) (3.47) −2S1(k)−k−1k(k+1)−ln(ξ)), B(0)3(ξ) = ∞∑k=0ξkΓ2(k+3)(k+1)Γ(2k+6)(2S1(k+2)−2S1(2k+5)−1k+1+ln(ξ)), (3.48) B(0)4(ξ) = −(ξ30+19)ζ2+(107ξ900+ζ26+1127)ln(ξ)−406713500ξ+136ln3(ξ) (3.49) −(ξ60+118)ln2(ξ)+ζ