This paper was converted on www.awesomepapers.org from LaTeX by an anonymous user.
Want to know more? Visit the Converter page.

aainstitutetext: School of Physics, Shandong University, Jinan, Shandong 250100, Chinabbinstitutetext: Center for High Energy Physics, Peking University, Beijing 100871, China

Application of the Meijer theorem in calculation of three-loop massive vacuum Feynman integrals and beyond

Jian Wang a    Dongyu Yang j.wang@sdu.edu.cn 202000141067@mail.sdu.edu.cn
Abstract

We present an analytical method to calculate the three-loop massive Feynman integral in arbitrary dimensions. The method is based on the Mellin-Barnes representation of the Feynman integral. The Meijer theorem and its corollary are used to perform the integration over the Gamma functions, exponential functions, and hypergeometric functions. We also discuss the application of the method in other multi-loop Feynman integrals.

1 Introduction

Multi-loop Feynman integrals play a crucial role in the application of quantum field theory. They are indispensable in calculating precise scattering cross-sections for collider processes that are important for testing the standard model. The massive vacuum loop integrals are one kind of the simplest multi-loop Feynman integrals. The results of these integrals are useful in understanding the structure of multi-loop integrals, e.g. the classes of constants and functions that would appear for Feynman loop integrals. They also find applications in the calculation of effective potentials of some theories Martin:2015eia . Recently, the method of differential equations Kotikov:1990kg ; Kotikov:1991pm has proven to be efficient in modern loop calculation and the vacuum integrals can be taken as the boundaries of the differential equations that more complex integrals satisfy.

The two-loop vacuum integrals with arbitrary masses have been obtained in full analytic form without expansion in the space-time parameter ϵ=(4−D)/2\epsilon=(4-D)/2 Davydychev:1992mt . The three-loop single-scale vacuum integrals have been calculated up to the finite part 𝒪​(ϵ0)\mathcal{O}(\epsilon^{0}) Broadhurst:1991fi ; Avdeev:1994db ; Broadhurst:1998rz ; Fleischer:1999mp ; Chetyrkin:1999qi , in which polylogarithms up to transcendental weight four are needed. The results were extended to weight six in Lee:2010hs ; Schroder:2005va ; Kniehl:2017ikj . The analytical results of three-loop two-scale vacuum integrals were available up to 𝒪​(ϵ0)\mathcal{O}(\epsilon^{0}) Bekavac:2009gz ; Grigo:2012ji . A special three-loop two-scale vacuum integral was calculated with full dependence on ϵ\epsilon Davydychev:2003mv . The finite parts of three-loop vacuum integrals with arbitrary mass pattern were computed numerically by using one- or two-dimensional integrals of elementary functions Freitas:2016zmy or by solving the coupled first-order differential equations Martin:2016bgz . Numerical results for even higher loop vacuum bubble diagrams of a specific type have also been presented Groote:2005ay . Recently, the analytical results have been investigated with the Gel’fand-Kapranov-Zelevinsky hypergeometric systems Gu:2018aya ; Gu:2020ypr ; Zhang:2023fil ; Zhang:2024mxd .

In this work, we present an analytical calculation of a three-loop massive Feynman integral in arbitrary dimensions. We adopted the Mellin-Barnes representation for the integral of Feynman parameters, and applied the Meijer theorem and its corollary to perform the integration over the Gamma functions, exponential functions, and hypergeometric functions. Our method is different from existing ones and we get a compact result of the three-loop vacuum integral. We also discuss the application of our method in other Feynman integrals.

2 Analytical calculation of the three-loop vacuum banana diagram

In this section, we use the well-known three-loop vacuum banana diagram with four massive propagators, labeled by BN, to demonstrate our approach. In subsection 2.1 we briefly review the definition of the Feynman integral of a three-loop vacuum diagram in DD dimensions and illustrate its Feynman parameter representation. In the following subsection 2.2 we introduce the Mellin-Barnes transformation for the Symanzik polynomials so that the integration over the Feynman parameters is easily performed. Then in subsection 2.3, we seek to lessen the number of parameters in Mellin-Barnes integrals and find it trivial to use Barnes’s first lemma to minimize a six-fold complex integral into a four-dimensional one. After that, we have to deal with integrals of a combination of different kinds of functions, including exponential functions, Γ\Gamma functions, and hypergeometric functions. To perform the integration, we propose an application of the Meijer theorem and its corollary.

2.1 Feynman representation and graph polynomials

k1k_{1}k2k_{2}k3k_{3}−k1+k2+k3-k_{1}+k_{2}+k_{3}
Figure 1: The three-loop massive vacuum diagram BN.

We begin our computation by writing down the normalized Feynman integral of the graph BN (see figure  1):

I≡e3​ϵ​γ​∫[d​k1]​[d​k2]​[d​k3]​1D1ν1​D2ν2​D3ν3​D4ν4I\equiv e^{3\epsilon\gamma}\int[\mathrm{d}k_{1}][\mathrm{d}k_{2}][\mathrm{d}k_{3}]\frac{1}{D_{1}^{\nu_{1}}D_{2}^{\nu_{2}}D_{3}^{\nu_{3}}D_{4}^{\nu_{4}}} (1)

with the Euler-Mascheroni constant γ=0.5772​⋯\gamma=0.5772\cdots, ν≡ν1+ν2+ν3+ν4\nu\equiv\nu_{1}+\nu_{2}+\nu_{3}+\nu_{4} and

[d​ki]\displaystyle[\mathrm{d}k_{i}] ≡dD​kii​πD/2,i=1,2,3,\displaystyle\equiv\frac{\mathrm{d}^{D}k_{i}}{i\pi^{D/2}},\qquad i=1,2,3,
Di\displaystyle D_{i} ≡−ki2+mi2,i=1,2,3,\displaystyle\equiv-k_{i}^{2}+m_{i}^{2},\qquad i=1,2,3,
D4\displaystyle D_{4} ≡−(−k1+k2+k3)2+m42.\displaystyle\equiv-(-k_{1}+k_{2}+k_{3})^{2}+m_{4}^{2}\,. (2)

The corresponding Symanzik polynomials 𝒰\mathcal{U} and ℱ\mathcal{F} are given by BOGNER2_2010

𝒰=a1​a2​a3+a1​a3​a4+a2​a3​a4+a1​a2​a4\mathcal{U}=a_{1}a_{2}a_{3}+a_{1}a_{3}a_{4}+a_{2}a_{3}a_{4}+a_{1}a_{2}a_{4} (3)

and

ℱ=(a1​a2​a3+a1​a3​a4+a2​a3​a4+a1​a2​a4)​(a1​m12+a2​m22+a3​m32+a4​m42).\mathcal{F}={(a_{1}a_{2}a_{3}+a_{1}a_{3}a_{4}+a_{2}a_{3}a_{4}+a_{1}a_{2}a_{4})(a_{1}m_{1}^{2}+a_{2}m_{2}^{2}+a_{3}m_{3}^{2}+a_{4}m_{4}^{2})}. (4)

The above integral is transformed to the integration over the Feynman parameters,

I=e3​ϵ​γ​Γ​(ν−3​D2)Γ​(ν1)​Γ​(ν2)​Γ​(ν3)​Γ​(ν4)​∫01da1​da2​da3​da4​δ​(1−a1−a2−a3−a4)×a1ν1−1​a2ν2−1​a3ν3−1​a4ν4−1​(a1​a2​a3+a1​a3​a4+a2​a3​a4+a1​a2​a4)−D2(a1​m12+a2​m22+a3​m32+a4​m42)ν−3​D2.I=\frac{e^{3\epsilon\gamma}\Gamma(\nu-\frac{3D}{2})}{\Gamma(\nu_{1})\Gamma(\nu_{2})\Gamma(\nu_{3})\Gamma(\nu_{4})}\int_{0}^{1}\mathrm{d}a_{1}\mathrm{d}a_{2}\mathrm{d}a_{3}\mathrm{d}a_{4}~{}\delta(1-a_{1}-a_{2}-a_{3}-a_{4})\\ \times a_{1}^{\nu_{1}-1}a_{2}^{\nu_{2}-1}a_{3}^{\nu_{3}-1}a_{4}^{\nu_{4}-1}\frac{(a_{1}a_{2}a_{3}+a_{1}a_{3}a_{4}+a_{2}a_{3}a_{4}+a_{1}a_{2}a_{4})^{-\frac{D}{2}}}{(a_{1}m_{1}^{2}+a_{2}m_{2}^{2}+a_{3}m_{3}^{2}+a_{4}m_{4}^{2})^{\nu-\frac{3D}{2}}}. (5)

2.2 Mellin-Barnes representation

Making use of the Mellin-Barnes transformation dubovyk2022mellinbarnes

(A1+A2+⋯+An)−c=1Γ​(c)​1(2​π​i)n−1​∫−i​∞i​∞dσ1​⋯​∫−i​∞i​∞dσn−1Γ​(−σ1)​⋯​Γ​(−σn−1)​Γ​(σ1+⋯+σn−1+c)​A1σ1​⋯​An−1σn−1​An−σ1−⋯​σn−1−c,(A_{1}+A_{2}+\cdots+A_{n})^{-c}=\frac{1}{\Gamma(c)}\frac{1}{(2\pi i)^{n-1}}\int_{-i\infty}^{i\infty}\mathrm{d}\sigma_{1}\cdots\int_{-i\infty}^{i\infty}\mathrm{d}\sigma_{n-1}\\ \Gamma(-\sigma_{1})\cdots\Gamma(-\sigma_{n-1})\Gamma(\sigma_{1}+\cdots+\sigma_{n-1}+c)A_{1}^{\sigma_{1}}\cdots A_{n-1}^{\sigma_{n-1}}A_{n}^{-\sigma_{1}-\cdots\sigma_{n-1}-c}, (6)

we derive

I=e3​ϵ​γΓ​(ν1)​Γ​(ν2)​Γ​(ν3)​Γ​(ν4)​Γ​(D/2)1(2​π​i)6×∫−i​∞i​∞dσ1​∫−i​∞i​∞dσ2​∫−i​∞i​∞dσ3​∫−i​∞i​∞dσ4​∫−i​∞i​∞dσ5​∫−i​∞i​∞dσ6Γ(−σ1)Γ(−σ2)Γ(−σ3)Γ(σ1+σ2+σ3+D2)Γ(−σ4)Γ(−σ5)Γ(−σ6)Γ(σ4+σ5+σ6+ν−3​D2)×(m12)σ4(m22)σ5(m32)σ6(m42)−σ4−σ5−σ6−ν+3​D2∫01da1da2da3da4⋅δ(1−a1−a2−a3−a4)×a1σ1+σ2+σ3+σ4+ν1−1​a2−σ3+σ5+ν2−D2−1​a3−σ2+σ6+ν3−D2−1​a4−σ1−σ4−σ5−σ6+ν4−ν+D−1.I=\frac{e^{3\epsilon\gamma}}{\Gamma(\nu_{1})\Gamma(\nu_{2})\Gamma(\nu_{3})\Gamma(\nu_{4})\Gamma(D/2)}\frac{1}{(2\pi i)^{6}}\times\\ \int_{-i\infty}^{i\infty}\mathrm{d}\sigma_{1}\int_{-i\infty}^{i\infty}\mathrm{d}\sigma_{2}\int_{-i\infty}^{i\infty}\mathrm{d}\sigma_{3}\int_{-i\infty}^{i\infty}\mathrm{d}\sigma_{4}\int_{-i\infty}^{i\infty}\mathrm{d}\sigma_{5}\int_{-i\infty}^{i\infty}\mathrm{d}\sigma_{6}\\ \Gamma(-\sigma_{1})\Gamma(-\sigma_{2})\Gamma(-\sigma_{3})\Gamma(\sigma_{1}+\sigma_{2}+\sigma_{3}+\frac{D}{2})\Gamma(-\sigma_{4})\Gamma(-\sigma_{5})\Gamma(-\sigma_{6})\Gamma(\sigma_{4}+\sigma_{5}+\sigma_{6}+\nu-\frac{3D}{2})\times\\ (m_{1}^{2})^{\sigma_{4}}(m_{2}^{2})^{\sigma_{5}}(m_{3}^{2})^{\sigma_{6}}(m_{4}^{2})^{-\sigma_{4}-\sigma_{5}-\sigma_{6}-\nu+\frac{3D}{2}}\int_{0}^{1}\mathrm{d}a_{1}\mathrm{d}a_{2}\mathrm{d}a_{3}\mathrm{d}a_{4}\cdot\delta(1-a_{1}-a_{2}-a_{3}-a_{4})\times\\ a_{1}^{\sigma_{1}+\sigma_{2}+\sigma_{3}+\sigma_{4}+\nu_{1}-1}a_{2}^{-\sigma_{3}+\sigma_{5}+\nu_{2}-\frac{D}{2}-1}a_{3}^{-\sigma_{2}+\sigma_{6}+\nu_{3}-\frac{D}{2}-1}a_{4}^{-\sigma_{1}-\sigma_{4}-\sigma_{5}-\sigma_{6}+\nu_{4}-\nu+D-1}. (7)

The integration over Feynman parameters can be done using the relation dubovyk2022mellinbarnes

∫01(∏j=1nd​aj⋅ajνj−1)​δ​(1−∑j=1naj)=∏j=1nΓ​(νj)Γ​(ν1+⋯+νn),\int_{0}^{1}\bigg{(}\prod_{j=1}^{n}\mathrm{d}a_{j}\cdot a_{j}^{\nu_{j}-1}\bigg{)}\delta\bigg{(}1-\sum_{j=1}^{n}a_{j}\bigg{)}=\frac{\prod_{j=1}^{n}\Gamma(\nu_{j})}{\Gamma(\nu_{1}+\cdots+\nu_{n})}, (8)

and we are left with

I=\displaystyle I= e3​ϵ​γ​(m42)−ν+3​D2Γ​(ν1)​Γ​(ν2)​Γ​(ν3)​Γ​(ν4)​Γ​(D/2)1(2​π​i)6×\displaystyle\frac{e^{3\epsilon\gamma}(m_{4}^{2})^{-\nu+\frac{3D}{2}}}{\Gamma(\nu_{1})\Gamma(\nu_{2})\Gamma(\nu_{3})\Gamma(\nu_{4})\Gamma(D/2)}\frac{1}{(2\pi i)^{6}}\times
∫−i​∞i​∞dσ1​∫−i​∞i​∞dσ2​∫−i​∞i​∞dσ3​∫−i​∞i​∞dσ4​∫−i​∞i​∞dσ5​∫−i​∞i​∞dσ6\displaystyle\int_{-i\infty}^{i\infty}\mathrm{d}\sigma_{1}\int_{-i\infty}^{i\infty}\mathrm{d}\sigma_{2}\int_{-i\infty}^{i\infty}\mathrm{d}\sigma_{3}\int_{-i\infty}^{i\infty}\mathrm{d}\sigma_{4}\int_{-i\infty}^{i\infty}\mathrm{d}\sigma_{5}\int_{-i\infty}^{i\infty}\mathrm{d}\sigma_{6}
Γ(−σ1)Γ(−σ2)Γ(−σ3)Γ(σ1+σ2+σ3+D2)Γ(−σ4)Γ(−σ5)Γ(−σ6)×\displaystyle\Gamma(-\sigma_{1})\Gamma(-\sigma_{2})\Gamma(-\sigma_{3})\Gamma(\sigma_{1}+\sigma_{2}+\sigma_{3}+\frac{D}{2})\Gamma(-\sigma_{4})\Gamma(-\sigma_{5})\Gamma(-\sigma_{6})\times
Γ(σ4+σ5+σ6+ν−3​D2)Γ(σ1+σ2+σ3+σ4+ν1)Γ(−σ3+σ5+ν2−D2)×\displaystyle\Gamma(\sigma_{4}+\sigma_{5}+\sigma_{6}+\nu-\frac{3D}{2})\Gamma(\sigma_{1}+\sigma_{2}+\sigma_{3}+\sigma_{4}+\nu_{1})\Gamma(-\sigma_{3}+\sigma_{5}+\nu_{2}-\frac{D}{2})\times
Γ(−σ2+σ6+ν3−D2)Γ(−σ1−σ4−σ5−σ6+ν4−ν+D)×\displaystyle\Gamma(-\sigma_{2}+\sigma_{6}+\nu_{3}-\frac{D}{2})\Gamma(-\sigma_{1}-\sigma_{4}-\sigma_{5}-\sigma_{6}+\nu_{4}-\nu+D)\times
(m12m42)σ4​(m22m42)σ5​(m32m42)σ6.\displaystyle\bigg{(}\frac{m_{1}^{2}}{m_{4}^{2}}\bigg{)}^{\sigma_{4}}\bigg{(}\frac{m_{2}^{2}}{m_{4}^{2}}\bigg{)}^{\sigma_{5}}\bigg{(}\frac{m_{3}^{2}}{m_{4}^{2}}\bigg{)}^{\sigma_{6}}. (9)

The integration contours must be chosen such that the poles of Γ​(⋯−σ)\Gamma(\cdots-\sigma) (Γ​(⋯+σ)\Gamma(\cdots+\sigma)) are to the right (left). For simplicity, we will focus on the case m1=m2=m3=m4=1m_{1}=m_{2}=m_{3}=m_{4}=1 in the following calculation, and discuss the general cases in the next section.

2.3 Reduction of integration parameters

In order to reduce the number of integration parameters, we apply Barnes’s first lemma barnes1908new .

Lemma 2.1 (Barnes’s first lemma).
12​π​i​∫−i​∞i​∞dσ​Γ​(a+σ)​Γ​(b+σ)​Γ​(c−σ)​Γ​(d−σ)=Γ​(a+c)​Γ​(a+d)​Γ​(b+c)​Γ​(b+d)Γ​(a+b+c+d).\frac{1}{2\pi i}\int_{-i\infty}^{i\infty}\mathrm{d}\sigma\Gamma(a+\sigma)\Gamma(b+\sigma)\Gamma(c-\sigma)\Gamma(d-\sigma)=\frac{\Gamma(a+c)\Gamma(a+d)\Gamma(b+c)\Gamma(b+d)}{\Gamma(a+b+c+d)}. (10)

The integral is convergent when Re​(a+b+c+d)<1\textrm{Re}(a+b+c+d)<1. However, this restriction can be removed after analytic continuation.

After integration over σ1\sigma_{1}, σ4\sigma_{4}, σ5\sigma_{5} and σ6\sigma_{6} in Eq. (9), we obtain

I=\displaystyle I= e3​ϵ​γΓ​(ν1)​Γ​(ν2)​Γ​(ν3)​Γ​(ν4)​Γ​(D/2)​1(2​π​i)2​∫−i​∞i​∞dσ2​dσ3\displaystyle\frac{e^{3\epsilon\gamma}}{\Gamma(\nu_{1})\Gamma(\nu_{2})\Gamma(\nu_{3})\Gamma(\nu_{4})\Gamma(D/2)}\frac{1}{(2\pi i)^{2}}\int_{-i\infty}^{i\infty}\mathrm{d}\sigma_{2}\mathrm{d}\sigma_{3}
Γ(−σ2)Γ(−σ3)Γ(σ2+σ3+D2)Γ(σ2+σ3+ν1)Γ(σ2+σ3+ν4)×\displaystyle\Gamma(-\sigma_{2})\Gamma(-\sigma_{3})\Gamma(\sigma_{2}+\sigma_{3}+\frac{D}{2})\Gamma(\sigma_{2}+\sigma_{3}+\nu_{1})\Gamma(\sigma_{2}+\sigma_{3}+\nu_{4})\times
Γ​(σ2+σ3−D2+ν1+ν4)​Γ​(−σ2−D2+ν3)​Γ​(−σ3−D2+ν2)Γ​(2​σ2+2​σ3+ν1+ν4).\displaystyle\frac{\Gamma(\sigma_{2}+\sigma_{3}-\frac{D}{2}+\nu_{1}+\nu_{4})\Gamma(-\sigma_{2}-\frac{D}{2}+\nu_{3})\Gamma(-\sigma_{3}-\frac{D}{2}+\nu_{2})}{\Gamma(2\sigma_{2}+2\sigma_{3}+\nu_{1}+\nu_{4})}. (11)

Without loss of generality, we set ν1=ν2=ν3=ν4=1\nu_{1}=\nu_{2}=\nu_{3}=\nu_{4}=1 below. Notice that the argument of the Gamma function in the denominator contains 2​σ22\sigma_{2} and 2​σ32\sigma_{3}. They can be normalized by using the Gauss multiplication formula abramowitz1968handbook .

Theorem 2.1 (Gauss multiplication formula).

For a positive integer nn,

(2​π)n−12​n12−n​z​Γ​(n​z)=∏k=0n−1Γ​(z+kn).(2\pi)^{\frac{n-1}{2}}n^{\frac{1}{2}-nz}\Gamma(nz)=\prod_{k=0}^{n-1}\Gamma\bigg{(}z+\frac{k}{n}\bigg{)}. (12)

Then we have

I=π​e3​ϵ​γ2​Γ​(2−ϵ)1(2​π​i)2∫−i​∞i​∞dσ2dσ3Γ(−σ2)Γ(−σ3)Γ(σ2+σ3+2−ϵ)Γ(σ2+σ3+1)×Γ​(σ2+σ3+ϵ)​Γ​(−σ2−1+ϵ)​Γ​(−σ3−1+ϵ)Γ​(σ2+σ3+32)⋅(14)σ2⋅(14)σ3.I=\frac{\sqrt{\pi}e^{3\epsilon\gamma}}{2\Gamma(2-\epsilon)}\frac{1}{(2\pi i)^{2}}\int_{-i\infty}^{i\infty}\mathrm{d}\sigma_{2}\mathrm{d}\sigma_{3}\Gamma(-\sigma_{2})\Gamma(-\sigma_{3})\Gamma(\sigma_{2}+\sigma_{3}+2-\epsilon)\Gamma(\sigma_{2}+\sigma_{3}+1)\times\\ \frac{\Gamma(\sigma_{2}+\sigma_{3}+\epsilon)\Gamma(-\sigma_{2}-1+\epsilon)\Gamma(-\sigma_{3}-1+\epsilon)}{\Gamma(\sigma_{2}+\sigma_{3}+\frac{3}{2})}\cdot\Big{(}\frac{1}{4}\Big{)}^{\sigma_{2}}\cdot\Big{(}\frac{1}{4}\Big{)}^{\sigma_{3}}. (13)

The integration over σ2\sigma_{2} or σ3\sigma_{3} can be considered as a general integral form of the hypergeometric function and can be carried out according to the Meijer theorem meijer1946 .

Theorem 2.2 (Meijer Theorem).

If the complex integral over s has a mixed product form of Γ\Gamma functions and exponential functions as the following

ℐ​(z)=12​π​i​∫−i​∞i​∞Γ​[(a)+s,(b)−s(c)+s,(d)−s]​zs​ds,\mathcal{I}(z)=\frac{1}{2\pi i}\int_{-i\infty}^{i\infty}\Gamma\Bigg{[}\begin{matrix}(a)+s,&(b)-s\\ (c)+s,&(d)-s\end{matrix}\Bigg{]}z^{s}\mathrm{d}s, (14)

which is the Meijer G-function, then the result is given by

ΣA(z)=∑μ=1Az−aμΓ[(a)μ′−aμ,(b)+aμ(c)−aμ,(d)+aμ]×FA+D−1B+C​[(b)+aμ,1+aμ−(c);1+aμ−(a)μ′,(d)+aμ;​(−1)A+C​z−1]\Sigma_{A}(z)=\sum_{\mu=1}^{A}z^{-a_{\mu}}\Gamma\Bigg{[}\begin{matrix}(a)_{\mu}^{\prime}-a_{\mu},&(b)+a_{\mu}\\ (c)-a_{\mu},&(d)+a_{\mu}\end{matrix}\Bigg{]}\times\\ {}_{B+C}F_{A+D-1}\Bigg{[}\begin{matrix}(b)+a_{\mu},1+a_{\mu}-(c);\\ 1+a_{\mu}-(a)_{\mu}^{\prime},(d)+a_{\mu};\end{matrix}(-1)^{A+C}z^{-1}\Bigg{]} (15)

or

ΣB(z)=∑ν=1BzbνΓ[(a)+bν,(b)ν′−bν(c)+bν,(d)−bν]×FB+C−1A+D​[(a)+bν,1+bν−(d);1+bν−(b)ν′,(c)+bν;​(−1)B+D​z],\Sigma_{B}(z)=\sum_{\nu=1}^{B}z^{b_{\nu}}\Gamma\Bigg{[}\begin{matrix}(a)+b_{\nu},&(b)_{\nu}^{\prime}-b_{\nu}\\ (c)+b_{\nu},&(d)-b_{\nu}\end{matrix}\Bigg{]}\times\\ {}_{A+D}F_{B+C-1}\Bigg{[}\begin{matrix}(a)+b_{\nu},1+b_{\nu}-(d);\\ 1+b_{\nu}-(b)_{\nu}^{\prime},(c)+b_{\nu};\end{matrix}(-1)^{B+D}z\Bigg{]}, (16)

where

Γ​[(a)+bν,(b)ν′−bν(c)+bν,(d)−bν]≡Γ​((a)+bν)​Γ​((b)ν′−bν)Γ​((c)+bν)​Γ​((d)−bν)\Gamma\Bigg{[}\begin{matrix}(a)+b_{\nu},&(b)_{\nu}^{\prime}-b_{\nu}\\ (c)+b_{\nu},&(d)-b_{\nu}\end{matrix}\Bigg{]}\equiv\frac{\Gamma((a)+b_{\nu})\Gamma((b)_{\nu}^{\prime}-b_{\nu})}{\Gamma((c)+b_{\nu})\Gamma((d)-b_{\nu})} (17)

and (a) is defined by a set of constants or parameters independent of s, i.e.,

(a)≡{a1,a2,…,aA}.(a)\equiv\{a_{1},a_{2},\dots,a_{A}\}. (18)

Meanwhile, the notation (a)ν′(a)_{\nu}^{\prime} represents a subset in which aνa_{\nu} is excluded, i.e.

(a)ν′≡{a1,…,aν−1,aν+1,…,aA},(a)_{\nu}^{\prime}\equiv\{a_{1},\dots,a_{\nu-1},a_{\nu+1},\dots,a_{A}\}, (19)

and A∈ℕA\in\mathbb{N} is the number of elements included in (a).

The choice between ΣA​(z)\Sigma_{A}(z) and ΣB​(z)\Sigma_{B}(z) is determined by a branch selection criterion:

  1. (i)

    ℐ​(z)=ΣA​(z)\mathcal{I}(z)=\Sigma_{A}(z)\quad when B+C<A+D\quad B+C<A+D  and  12​π​|A+B−C−D|>|arg⁡z|\frac{1}{2}\pi|A+B-C-D|>|\arg z|

    or B+C=A+D\quad B+C=A+D\quad and |z|>1\quad|z|>1.

  2. (ii)

    ℐ​(z)=ΣB​(z)\mathcal{I}(z)=\Sigma_{B}(z)\quad when B+C>A+D\quad B+C>A+D  and  12​π​|A+B−C−D|>|arg⁡z|\frac{1}{2}\pi|A+B-C-D|>|\arg z|

    or B+C=A+D\quad B+C=A+D\quad and |z|<1\quad|z|<1.

  3. (iii)

    ℐ​(1)=ΣA​(1)=ΣB​(1)\mathcal{I}(1)=\Sigma_{A}(1)=\Sigma_{B}(1)\quad when A−C=B−D≥0\quad A-C=B-D\geq 0 and Re​(Σ​c+Σ​d−Σ​a−Σ​b)>0\textrm{Re}(\Sigma c+\Sigma d-\Sigma a-\Sigma b)>0.

After applying the Meijer theorem, eq. (13) is converted to

I=π​e3​ϵ​γ​Γ​(−1+ϵ)2​Γ​(2−ϵ)12​π​i∫−i​∞i​∞dσ2⋅Γ[σ2+2−ϵ,σ2+1,σ2+ϵ,−σ2,−σ2−1+ϵσ2+32]×F23​[σ2+2−ϵ,σ2+1,σ2+ϵ;2−ϵ,σ2+32;14]​(14)σ2+π​e3​ϵ​γ​Γ​(1−ϵ)2​Γ​(2−ϵ)12​π​i⋅(14)−1+ϵ∫−i​∞i​∞dσ2⋅Γ[σ2+1,σ2+ϵ,σ2+2​ϵ−1,−σ2,−σ2−1+ϵσ2+12+ϵ]×F23​[σ2−1+2​ϵ,σ2+1,σ2+ϵ;ϵ,σ2+12+ϵ;14]​(14)σ2.I=\frac{\sqrt{\pi}e^{3\epsilon\gamma}\Gamma(-1+\epsilon)}{2\Gamma(2-\epsilon)}\frac{1}{2\pi i}\int_{-i\infty}^{i\infty}\mathrm{d}\sigma_{2}\cdot\Gamma\Bigg{[}\begin{matrix}\sigma_{2}+2-\epsilon,\sigma_{2}+1,\sigma_{2}+\epsilon,-\sigma_{2},-\sigma_{2}-1+\epsilon\\ \sigma_{2}+\dfrac{3}{2}\end{matrix}\Bigg{]}\times\\ {}_{3}F_{2}\Bigg{[}\begin{matrix}\sigma_{2}+2-\epsilon,\sigma_{2}+1,\sigma_{2}+\epsilon;\\ 2-\epsilon,\sigma_{2}+\dfrac{3}{2};\end{matrix}\quad\frac{1}{4}\Bigg{]}\bigg{(}\frac{1}{4}\bigg{)}^{\sigma_{2}}+\\ \frac{\sqrt{\pi}e^{3\epsilon\gamma}\Gamma(1-\epsilon)}{2\Gamma(2-\epsilon)}\frac{1}{2\pi i}\cdot\bigg{(}\frac{1}{4}\bigg{)}^{-1+\epsilon}\int_{-i\infty}^{i\infty}\mathrm{d}\sigma_{2}\cdot\Gamma\Bigg{[}\begin{matrix}\sigma_{2}+1,\sigma_{2}+\epsilon,\sigma_{2}+2\epsilon-1,-\sigma_{2},-\sigma_{2}-1+\epsilon\\ \sigma_{2}+\dfrac{1}{2}+\epsilon\end{matrix}\Bigg{]}\times\\ {}_{3}F_{2}\Bigg{[}\begin{matrix}\sigma_{2}-1+2\epsilon,\sigma_{2}+1,\sigma_{2}+\epsilon;\\ \epsilon,\sigma_{2}+\dfrac{1}{2}+\epsilon;\end{matrix}\quad\frac{1}{4}\Bigg{]}\bigg{(}\frac{1}{4}\bigg{)}^{\sigma_{2}}. (20)

The integrands of eq. (20) appear to be a compound of Γ\Gamma functions, hypergeometric functions and exponential functions. From the Meijer theorem, we derive the following corollary which helps to perform the integration111A similar form of this corollary was presented by L. J. Slater in  slater1966generalized . However, the result contains some typos.. The details of the proof are provided in appendix A.

Corollary 2.1.

If the complex integral over s has a mixed product form of Γ\Gamma functions, exponential functions, and hypergeometric functions as the following

ℐ(z)=12​π​i∫−i​∞i​∞Γ[(a)+s,(b)−s,(g)+s,(h)−s(c)+s,(d)−s,(j)+s,(k)−s]×FC+D+FA+B+E​[(a)+s,(b)−s,(e);(c)+s,(d)−s,(f);​x]​zs​d​s,\mathcal{I}(z)=\frac{1}{2\pi i}\int_{-i\infty}^{i\infty}\Gamma\Bigg{[}\begin{matrix}(a)+s,&(b)-s,&(g)+s,(h)-s\\ (c)+s,&(d)-s,&(j)+s,(k)-s\end{matrix}\Bigg{]}\times\\ {}_{A+B+E}F_{C+D+F}\Bigg{[}\begin{matrix}(a)+s,&(b)-s,&(e);\\ (c)+s,&(d)-s,&(f);\end{matrix}~{}x\Bigg{]}z^{s}\mathrm{d}s, (21)

then the result is given by

ΣA(z)=∑μ=1AΓ[(a)μ′−aμ,(b)+aμ,(g)−aμ,(h)+aμ(c)−aμ,(d)+aμ,(j)−aμ,(k)+aμ]×∑m=0∞∑n=0∞((b)+aμ)2​m+n​((h)+aμ)m+n​((e))m​(1+aμ−(c))n​(1+aμ−(j))m+n(1+aμ−(a)μ′)n​(1+aμ−(g))m+n​((f))m​((d)+aμ)2​m+n​((k)+aμ)m+n⋅xm​z−aμ−m−n​(−1)n​(A+G−C−J)+m​(G−J)m!​n!+∑μ=1GΓ[(a)−gμ,(b)+gμ,(g)μ′−gμ,(h)+gμ(c)−gμ,(d)+gμ,(j)−gμ,(k)+gμ]×∑m=0∞∑n=0∞((a)−gμ)m−n​((b)+gμ)m+n​(1+gμ−(j))n​((h)+gμ)n​((e))m((c)−gμ)m−n​((d)+gμ)m+n​(1+gμ−(g)μ′)n​((k)+gμ)n​((f))m⋅xm​z−gμ−n​(−1)n​(G−J)m!​n!,\Sigma_{A}(z)=\sum_{\mu=1}^{A}\Gamma\Bigg{[}\begin{matrix}(a)_{\mu}^{\prime}-a_{\mu},&(b)+a_{\mu},&(g)-a_{\mu},&(h)+a_{\mu}\\ (c)-a_{\mu},&(d)+a_{\mu},&(j)-a_{\mu},&(k)+a_{\mu}\end{matrix}\Bigg{]}\times\\ \sum_{m=0}^{\infty}\sum_{n=0}^{\infty}\frac{((b)+a_{\mu})_{2m+n}((h)+a_{\mu})_{m+n}((e))_{m}(1+a_{\mu}-(c))_{n}(1+a_{\mu}-(j))_{m+n}}{(1+a_{\mu}-(a)_{\mu}^{\prime})_{n}(1+a_{\mu}-(g))_{m+n}((f))_{m}((d)+a_{\mu})_{2m+n}((k)+a_{\mu})_{m+n}}\\ \cdot\frac{x^{m}z^{-a_{\mu}-m-n}(-1)^{n(A+G-C-J)+m(G-J)}}{m!n!}\\ +\sum_{\mu=1}^{G}\Gamma\Bigg{[}\begin{matrix}(a)-g_{\mu},&(b)+g_{\mu},&(g)_{\mu}^{\prime}-g_{\mu},&(h)+g_{\mu}\\ (c)-g_{\mu},&(d)+g_{\mu},&(j)-g_{\mu},&(k)+g_{\mu}\end{matrix}\Bigg{]}\times\\ \sum_{m=0}^{\infty}\sum_{n=0}^{\infty}\frac{((a)-g_{\mu})_{m-n}((b)+g_{\mu})_{m+n}(1+g_{\mu}-(j))_{n}((h)+g_{\mu})_{n}((e))_{m}}{((c)-g_{\mu})_{m-n}((d)+g_{\mu})_{m+n}(1+g_{\mu}-(g)_{\mu}^{\prime})_{n}((k)+g_{\mu})_{n}((f))_{m}}\\ \cdot\frac{x^{m}z^{-g_{\mu}-n}(-1)^{n(G-J)}}{m!n!}, (22)

or

ΣB(z)=∑ν=1BΓ[(a)+bν,(b)ν′−bν,(g)+bν,(h)−bν(c)+bν,(d)−bν,(j)+bν,(k)−bν]×∑m=0∞∑n=0∞((a)+bν)2​m+n​((g)+bν)m+n​((e))m​(1+bν−(d))n​(1+bν−(k))m+n(1+bν−(b)ν′)n​(1+bν−(h))m+n​((f))m​((c)+bν)2​m+n​((j)+bν)m+n⋅xm​zbν+m+n​(−1)n​(B+H−K−D)+m​(H−K)m!​n!+∑ν=1HΓ[(a)+hν,(b)−hν,(g)+hν,(h)ν′−hν(c)+hν,(d)−hν,(j)+hν,(k)−hν]×∑m=0∞∑n=0∞((a)+hν)m+n​((b)−hν)m−n​(1+hν−(k))n​((g)+hν)n​((e))m((d)−hν)m−n​((c)+hν)m+n​(1+hν−(h)ν′)n​((j)+hν)n​((f))m⋅xm​zhν+n​(−1)n​(H−K)m!​n!,\Sigma_{B}(z)=\sum_{\nu=1}^{B}\Gamma\Bigg{[}\begin{matrix}(a)+b_{\nu},&(b)_{\nu}^{\prime}-b_{\nu},&(g)+b_{\nu},&(h)-b_{\nu}\\ (c)+b_{\nu},&(d)-b_{\nu},&(j)+b_{\nu},&(k)-b_{\nu}\end{matrix}\Bigg{]}\times\\ \sum_{m=0}^{\infty}\sum_{n=0}^{\infty}\frac{((a)+b_{\nu})_{2m+n}((g)+b_{\nu})_{m+n}((e))_{m}(1+b_{\nu}-(d))_{n}(1+b_{\nu}-(k))_{m+n}}{(1+b_{\nu}-(b)_{\nu}^{\prime})_{n}(1+b_{\nu}-(h))_{m+n}((f))_{m}((c)+b_{\nu})_{2m+n}((j)+b_{\nu})_{m+n}}\\ \cdot\frac{x^{m}z^{b_{\nu}+m+n}(-1)^{n(B+H-K-D)+m(H-K)}}{m!n!}\\ +\sum_{\nu=1}^{H}\Gamma\Bigg{[}\begin{matrix}(a)+h_{\nu},&(b)-h_{\nu},&(g)+h_{\nu},&(h)_{\nu}^{\prime}-h_{\nu}\\ (c)+h_{\nu},&(d)-h_{\nu},&(j)+h_{\nu},&(k)-h_{\nu}\end{matrix}\Bigg{]}\times\\ \sum_{m=0}^{\infty}\sum_{n=0}^{\infty}\frac{((a)+h_{\nu})_{m+n}((b)-h_{\nu})_{m-n}(1+h_{\nu}-(k))_{n}((g)+h_{\nu})_{n}((e))_{m}}{((d)-h_{\nu})_{m-n}((c)+h_{\nu})_{m+n}(1+h_{\nu}-(h)_{\nu}^{\prime})_{n}((j)+h_{\nu})_{n}((f))_{m}}\\ \cdot\frac{x^{m}z^{h_{\nu}+n}(-1)^{n(H-K)}}{m!n!}, (23)

where the Pochhammer symbol (a)n(a)_{n} knuth1992notes ; abramowitz1968handbook is defined by

(a)n≡Γ​(a+n)Γ​(a),(a)_{n}\equiv\frac{\Gamma(a+n)}{\Gamma(a)}, (24)

and the hypergeometric function FC+D+FA+B+E​(x){}_{A+B+E}F_{C+D+F}(x) is absolutely and uniformly convergent in xx.

Provided that 12​π​|A+G+B+H−C−D−J−K|>|arg⁡z|\frac{1}{2}\pi|A+G+B+H-C-D-J-K|>|\arg z|, the choice between ΣA​(z)\Sigma_{A}(z) and ΣB​(z)\Sigma_{B}(z) is determined by a branch selection criterion:

  1. (i)

    ℐ​(z)=ΣA​(z)\mathcal{I}(z)=\Sigma_{A}(z)   when A+G+D+K>B+H+C+J\qquad A+G+D+K>B+H+C+J\qquad

    or A+G+D+K=B+H+C+J\qquad A+G+D+K=B+H+C+J   and   |z|>1|z|>1.

  2. (ii)

    ℐ​(z)=ΣB​(z)\mathcal{I}(z)=\Sigma_{B}(z)   when A+G+D+K<B+H+C+J\qquad A+G+D+K<B+H+C+J\qquad

    or A+G+D+K=B+H+C+J\qquad A+G+D+K=B+H+C+J   and   |z|<1|z|<1.

Also, provided that z=1z=1, and Re(Σc+Σd+Σj+Σk−+Σa−Σb−Σg−Σh)>0\textrm{Re}(\Sigma c+\Sigma d+\Sigma j+\Sigma k-+\Sigma a-\Sigma b-\Sigma g-\Sigma h)>0,

  1. (iii)

    ℐ​(1)=ΣA​(1)=ΣB​(1)\mathcal{I}(1)=\Sigma_{A}(1)=\Sigma_{B}(1)  when A+G−C−J=B+H−D−K≥0\qquad A+G-C-J=B+H-D-K\geq 0.

Finally, we obtain the analytical result for the three-loop massive vacuum integral BN in arbitrary dimensions:

I=π​e3​γ​ϵ​Γ2​(−1+ϵ)​Γ​(ϵ)2​Γ​(32)​∑m=0∞∑n=0∞(2−ϵ)m+n​(1)m+n​(ϵ)m+n(2−ϵ)n​(32)m+n​(2−ϵ)m⋅1m!​n!⋅(14)m+n+π​e3​γ​ϵ​Γ​(ϵ)​Γ​(−1+2​ϵ)​Γ​(1−ϵ)​Γ​(−1+ϵ)Γ​(2−ϵ)​Γ​(12+ϵ)⋅(14)−1+ϵ×∑m=0∞∑n=0∞(1)m+n​(ϵ)m+n​(−1+2​ϵ)m+n(ϵ)n​(12+ϵ)m+n​(2−ϵ)m⋅1m!​n!⋅(14)m+n+π​e3​γ​ϵ​Γ​(ϵ)​Γ​(−1+2​ϵ)​Γ​(−2+3​ϵ)​Γ2​(1−ϵ)2​Γ​(2−ϵ)​Γ​(−12+2​ϵ)⋅(14)−2+2​ϵ×∑m=0∞∑n=0∞(ϵ)m+n​(−1+2​ϵ)m+n​(−2+3​ϵ)m+n(ϵ)n​(−12+2​ϵ)m+n​(ϵ)m⋅1m!​n!⋅(14)m+n.I=\frac{\sqrt{\pi}e^{3\gamma\epsilon}\Gamma^{2}(-1+\epsilon)\Gamma(\epsilon)}{2\Gamma(\frac{3}{2})}\sum_{m=0}^{\infty}\sum_{n=0}^{\infty}\frac{(2-\epsilon)_{m+n}(1)_{m+n}(\epsilon)_{m+n}}{(2-\epsilon)_{n}(\frac{3}{2})_{m+n}(2-\epsilon)_{m}}\cdot\frac{1}{m!n!}\cdot\bigg{(}\frac{1}{4}\bigg{)}^{m+n}+\\ \frac{\sqrt{\pi}e^{3\gamma\epsilon}\Gamma(\epsilon)\Gamma(-1+2\epsilon)\Gamma(1-\epsilon)\Gamma(-1+\epsilon)}{\Gamma(2-\epsilon)\Gamma(\frac{1}{2}+\epsilon)}\cdot\bigg{(}\frac{1}{4}\bigg{)}^{-1+\epsilon}\times\\ \sum_{m=0}^{\infty}\sum_{n=0}^{\infty}\frac{(1)_{m+n}(\epsilon)_{m+n}(-1+2\epsilon)_{m+n}}{(\epsilon)_{n}(\frac{1}{2}+\epsilon)_{m+n}(2-\epsilon)_{m}}\cdot\frac{1}{m!n!}\cdot\bigg{(}\frac{1}{4}\bigg{)}^{m+n}+\\ \frac{\sqrt{\pi}e^{3\gamma\epsilon}\Gamma(\epsilon)\Gamma(-1+2\epsilon)\Gamma(-2+3\epsilon)\Gamma^{2}(1-\epsilon)}{2\Gamma(2-\epsilon)\Gamma(-\frac{1}{2}+2\epsilon)}\cdot\bigg{(}\frac{1}{4}\bigg{)}^{-2+2\epsilon}\times\\ \sum_{m=0}^{\infty}\sum_{n=0}^{\infty}\frac{(\epsilon)_{m+n}(-1+2\epsilon)_{m+n}(-2+3\epsilon)_{m+n}}{(\epsilon)_{n}(-\frac{1}{2}+2\epsilon)_{m+n}(\epsilon)_{m}}\cdot\frac{1}{m!n!}\cdot\bigg{(}\frac{1}{4}\bigg{)}^{m+n}. (25)

This result agrees with that in ref. Gu:2018aya where a different method has been adopted. The series of expansion around ϵ=0\epsilon=0 is given by

I=2ϵ−3+223ϵ−2+(272+8ln2+π22−∂2∂x​∂yF12[x,y;−1/2;1]|x=−1/2,y=−2−2∂2∂x​∂yF12[−1/2,x;y;1]|x=−2,y=−1/2)ϵ−1+𝒪(ϵ0),I=2\epsilon^{-3}+\frac{22}{3}\epsilon^{-2}+\bigg{(}\frac{27}{2}+8\ln 2+\frac{\pi^{2}}{2}-\frac{\partial^{2}}{\partial x\partial y}{}_{2}F_{1}\Big{[}\begin{matrix}x,\ y;\\ -1/2;\end{matrix}\quad 1\Big{]}\Big{|}_{x=-1/2,y=-2}\\ -2\frac{\partial^{2}}{\partial x\partial y}{}_{2}F_{1}\Big{[}\begin{matrix}-1/2,\ x;\\ y;\end{matrix}\quad 1\Big{]}\Big{|}_{x=-2,y=-1/2}\bigg{)}\epsilon^{-1}+\mathcal{O}(\epsilon^{0}), (26)

which completely coincides with the numerical calculations in kniehl2017three .

3 Application of the method in other Feynman integrals

The application of our calculation method in the above example is successful. In this section, we present several more examples and a general algorithm for the calculation of Feynman integrals.

BN diagrams with arbitrary masses

In the above calculation, we have made an assumption that all masses of the four propagators are identical. This does not mean that our method is only applicable in such a simple case. Now we illustrate how to perform calculations for the BN diagrams with different masses.

Let us start from eq. (9). It is easy to integrate over σ1\sigma_{1} using Barnes’s first lemma, and we get

I=e3​ϵ​γ​(m42)−ν+3​D2Γ​(ν1)​Γ​(ν2)​Γ​(ν3)​Γ​(ν4)​Γ​(D/2)1(2​π​i)5×∫−i​∞i​∞dσ2∫−i​∞i​∞dσ3∫−i​∞i​∞dσ4∫−i​∞i​∞dσ5∫−i​∞i​∞dσ6(m12m42)σ4(m22m42)σ5(m32m42)σ6×Γ(−σ2)Γ(−σ3)Γ(−σ4)Γ(−σ5)Γ(−σ6)Γ(σ4+σ5+σ6+ν−3​D2)×Γ(−σ3+σ5+ν2−D2)Γ(−σ2+σ6+ν3−D2)Γ(σ2+σ3−σ4−σ5−σ6+3​D2+ν4−ν)×Γ​(σ2+σ3+D2)​Γ​(σ2+σ3+σ4+ν1)​Γ​(σ2+σ3−σ5−σ6+D+ν4+ν1−ν)Γ​(2​σ2+2​σ3−σ5−σ6+3​D2+ν1+ν4−ν).I=\frac{e^{3\epsilon\gamma}(m_{4}^{2})^{-\nu+\frac{3D}{2}}}{\Gamma(\nu_{1})\Gamma(\nu_{2})\Gamma(\nu_{3})\Gamma(\nu_{4})\Gamma(D/2)}\frac{1}{(2\pi i)^{5}}\times\\ \int_{-i\infty}^{i\infty}\mathrm{d}\sigma_{2}\int_{-i\infty}^{i\infty}\mathrm{d}\sigma_{3}\int_{-i\infty}^{i\infty}\mathrm{d}\sigma_{4}\int_{-i\infty}^{i\infty}\mathrm{d}\sigma_{5}\int_{-i\infty}^{i\infty}\mathrm{d}\sigma_{6}\bigg{(}\frac{m_{1}^{2}}{m_{4}^{2}}\bigg{)}^{\sigma_{4}}\bigg{(}\frac{m_{2}^{2}}{m_{4}^{2}}\bigg{)}^{\sigma_{5}}\bigg{(}\frac{m_{3}^{2}}{m_{4}^{2}}\bigg{)}^{\sigma_{6}}\times\\ \Gamma(-\sigma_{2})\Gamma(-\sigma_{3})\Gamma(-\sigma_{4})\Gamma(-\sigma_{5})\Gamma(-\sigma_{6})\Gamma(\sigma_{4}+\sigma_{5}+\sigma_{6}+\nu-\frac{3D}{2})\times\\ \Gamma(-\sigma_{3}+\sigma_{5}+\nu_{2}-\frac{D}{2})\Gamma(-\sigma_{2}+\sigma_{6}+\nu_{3}-\frac{D}{2})\Gamma(\sigma_{2}+\sigma_{3}-\sigma_{4}-\sigma_{5}-\sigma_{6}+\frac{3D}{2}+\nu_{4}-\nu)\times\\ \frac{\Gamma(\sigma_{2}+\sigma_{3}+\frac{D}{2})\Gamma(\sigma_{2}+\sigma_{3}+\sigma_{4}+\nu_{1})\Gamma(\sigma_{2}+\sigma_{3}-\sigma_{5}-\sigma_{6}+D+\nu_{4}+\nu_{1}-\nu)}{\Gamma(2\sigma_{2}+2\sigma_{3}-\sigma_{5}-\sigma_{6}+\frac{3D}{2}+\nu_{1}+\nu_{4}-\nu)}. (27)

Assuming

m1≥max⁡(m2,m3,m4),m_{1}\geq\max(m_{2},\ m_{3},\ m_{4}), (28)

the Meijer theorem can be applied to the integration of σ4\sigma_{4}, yielding

I=e3​ϵ​γ​(m42)−ν+3​D2Γ​(ν1)​Γ​(ν2)​Γ​(ν3)​Γ​(ν4)​Γ​(D/2)​1(2​π​i)4​∫−i​∞i​∞dσ2​dσ3​dσ5​dσ6×Γ​(−σ2)​Γ​(−σ3)​Γ​(−σ5)​Γ​(−σ6)​Γ​(−σ3+σ5+ν2−D/2)​Γ​(−σ2+σ6+ν3−D/2)×Γ​(σ2+σ3+D/2)​Γ​(σ2+σ3−σ5−σ6+D+ν4+ν1−ν)Γ​(2​σ2+2​σ3−σ5−σ6+3​D2+ν1+ν4−ν)​(m22m42)σ5​(m32m42)σ6×{(m12m42)−σ2−σ3−ν1Γ(−σ2−σ3+σ5+σ6+ν−ν1−3D/2)Γ(σ2+σ3+ν1)×Γ​(2​σ2+2​σ3−σ5−σ6+3​D/2+ν1+ν4−ν)×F12​[σ2+σ3+ν1,2​σ2+2​σ3−σ5−σ6+3​D2+ν1+ν4−ν;σ2+σ3−σ5−σ6−ν+ν1+3​D2+1;m42m12]+(m12m42)−σ5−σ6−ν+3​D2​Γ​(σ2+σ3−σ5−σ6−ν+ν1+3​D/2)​Γ​(σ5+σ6+ν−3​D/2)×Γ(σ2+σ3+ν4)2F1[σ5+σ6+ν−3​D2,σ2+σ3+ν4;−σ2−σ3+σ5+σ6+ν−ν1−3​D2+1;m42m12]}.I=\frac{e^{3\epsilon\gamma}(m_{4}^{2})^{-\nu+\frac{3D}{2}}}{\Gamma(\nu_{1})\Gamma(\nu_{2})\Gamma(\nu_{3})\Gamma(\nu_{4})\Gamma(D/2)}\frac{1}{(2\pi i)^{4}}\int_{-i\infty}^{i\infty}\mathrm{d}\sigma_{2}\mathrm{d}\sigma_{3}\mathrm{d}\sigma_{5}\mathrm{d}\sigma_{6}\\ \times\Gamma(-\sigma_{2})\Gamma(-\sigma_{3})\Gamma(-\sigma_{5})\Gamma(-\sigma_{6})\Gamma(-\sigma_{3}+\sigma_{5}+\nu_{2}-D/2)\Gamma(-\sigma_{2}+\sigma_{6}+\nu_{3}-D/2)\\ \times\frac{\Gamma(\sigma_{2}+\sigma_{3}+D/2)\Gamma(\sigma_{2}+\sigma_{3}-\sigma_{5}-\sigma_{6}+D+\nu_{4}+\nu_{1}-\nu)}{\Gamma(2\sigma_{2}+2\sigma_{3}-\sigma_{5}-\sigma_{6}+\frac{3D}{2}+\nu_{1}+\nu_{4}-\nu)}\bigg{(}\frac{m_{2}^{2}}{m_{4}^{2}}\bigg{)}^{\sigma_{5}}\bigg{(}\frac{m_{3}^{2}}{m_{4}^{2}}\bigg{)}^{\sigma_{6}}\\ \times\Bigg{\{}\bigg{(}\frac{m_{1}^{2}}{m_{4}^{2}}\bigg{)}^{-\sigma_{2}-\sigma_{3}-\nu_{1}}\Gamma(-\sigma_{2}-\sigma_{3}+\sigma_{5}+\sigma_{6}+\nu-\nu_{1}-3D/2)\Gamma(\sigma_{2}+\sigma_{3}+\nu_{1})\\ \times\Gamma(2\sigma_{2}+2\sigma_{3}-\sigma_{5}-\sigma_{6}+3D/2+\nu_{1}+\nu_{4}-\nu)\\ \times{}_{2}F_{1}\Bigg{[}\begin{matrix}\sigma_{2}+\sigma_{3}+\nu_{1},2\sigma_{2}+2\sigma_{3}-\sigma_{5}-\sigma_{6}+\frac{3D}{2}+\nu_{1}+\nu_{4}-\nu;\\ \sigma_{2}+\sigma_{3}-\sigma_{5}-\sigma_{6}-\nu+\nu_{1}+\frac{3D}{2}+1;\end{matrix}\quad\frac{m_{4}^{2}}{m_{1}^{2}}\Bigg{]}\\ +\bigg{(}\frac{m_{1}^{2}}{m_{4}^{2}}\bigg{)}^{-\sigma_{5}-\sigma_{6}-\nu+\frac{3D}{2}}\Gamma(\sigma_{2}+\sigma_{3}-\sigma_{5}-\sigma_{6}-\nu+\nu_{1}+3D/2)\Gamma(\sigma_{5}+\sigma_{6}+\nu-3D/2)\\ \times\Gamma(\sigma_{2}+\sigma_{3}+\nu_{4})_{2}F_{1}\Bigg{[}\begin{matrix}\sigma_{5}+\sigma_{6}+\nu-\frac{3D}{2},\ \sigma_{2}+\sigma_{3}+\nu_{4};\\ -\sigma_{2}-\sigma_{3}+\sigma_{5}+\sigma_{6}+\nu-\nu_{1}-\frac{3D}{2}+1;\end{matrix}\quad\frac{m_{4}^{2}}{m_{1}^{2}}\Bigg{]}\Bigg{\}}. (29)

It is manifest that our next step is to integrate over σ5\sigma_{5} (or σ6\sigma_{6}) using Corollary 2.1. The result would contain Pochhammer symbols which are not suitable for the application of the Meijer theorem again. Therefore it is important to rewrite ΣA\Sigma_{A} and ΣB\Sigma_{B} in Corollary 2.1 in terms of Γ\Gamma functions:

ΣA(z)=∑μ=1A∑m=0∞∑n=0∞Γ[(a)μ′−aμ,1+aμ−(a)μ′,(g)−aμ,1+aμ−(g),(f),(c)−aμ,1+aμ−(c),(j)−aμ,1−aμ−(j),(e),(b)+aμ+2​m+n,(h)+aμ+m+n,1+aμ−(c)+n,1+aμ−(j)+m+n,(e)+m(d)+aμ+2​m+n,(k)+aμ+m+n,1+aμ−(a)μ′+n,1+aμ−(g)+m+n,(f)+m]×xm​z−aμ−m−n​(−1)n​(A+G−C−J)+m​(G−J)m!​n!+∑μ=1G∑m=0∞∑n=0∞Γ[(g)μ′−gμ,1+gμ−(g)μ′,(f),(a)−gμ+m−n,(j)−gμ,1+gμ−(j),(e),(c)−gμ+m−n,(b)+gμ+m+n,1+gμ−(j)+n,(h)+gμ+n,(e)+m(d)+gμ+m+n,1+gμ−(g)μ′+n,(k)+gμ+n,(f)+m]xm​z−gμ−n​(−1)n​(G−J)m!​n!,\Sigma_{A}(z)=\sum_{\mu=1}^{A}\sum_{m=0}^{\infty}\sum_{n=0}^{\infty}\Gamma\Bigg{[}\begin{matrix}(a)_{\mu}^{\prime}-a_{\mu},&1+a_{\mu}-(a)_{\mu}^{\prime},&(g)-a_{\mu},&1+a_{\mu}-(g),&(f),\\ (c)-a_{\mu},&1+a_{\mu}-(c),&(j)-a_{\mu},&1-a_{\mu}-(j),&(e),\end{matrix}\\ \begin{matrix}(b)+a_{\mu}+2m+n,&(h)+a_{\mu}+m+n,&1+a_{\mu}-(c)+n,&1+a_{\mu}-(j)+m+n,&(e)+m\\ (d)+a_{\mu}+2m+n,&(k)+a_{\mu}+m+n,&1+a_{\mu}-(a)_{\mu}^{\prime}+n,&1+a_{\mu}-(g)+m+n,&(f)+m\end{matrix}\Bigg{]}\\ \times\frac{x^{m}z^{-a_{\mu}-m-n}(-1)^{n(A+G-C-J)+m(G-J)}}{m!n!}\\ +\sum_{\mu=1}^{G}\sum_{m=0}^{\infty}\sum_{n=0}^{\infty}\Gamma\Bigg{[}\begin{matrix}(g)_{\mu}^{\prime}-g_{\mu},&1+g_{\mu}-(g)_{\mu}^{\prime},&(f),&(a)-g_{\mu}+m-n,\\ (j)-g_{\mu},&1+g_{\mu}-(j),&(e),&(c)-g_{\mu}+m-n,\end{matrix}\\ \begin{matrix}(b)+g_{\mu}+m+n,&1+g_{\mu}-(j)+n,&(h)+g_{\mu}+n,&(e)+m\\ (d)+g_{\mu}+m+n,&1+g_{\mu}-(g)_{\mu}^{\prime}+n,&(k)+g_{\mu}+n,&(f)+m\end{matrix}\Bigg{]}\frac{x^{m}z^{-g_{\mu}-n}(-1)^{n(G-J)}}{m!n!}, (30)
ΣB(z)=∑ν=1B∑m=0∞∑n=0∞Γ[(b)ν′−bν,(h)−bν,1+bν−(b)ν′,1+bν−(h),(f),(d)−bν,(k)−bν,1+bν−(d),1+bν−(k),(e),(a)+bν+2​m+n,(g)+bν+m+n,(e)+m,1+bν−(d)+n,1+bν−(k)+m+n(c)+bν+2​m+n,(j)+bν+m+n,(f)+m,1+bν−(b)ν′+n,1+bν−(h)+m+n]×xm​zbν+m+n​(−1)n​(B+H−K−D)+m​(H−K)m!​n!+∑ν=1H∑m=0∞∑n=0∞Γ[(h)ν′−hν,1+hν−(h)ν′,(f),(a)+hν+m+n,(b)−hν+m−n,(k)−hν,1+hν−(k),(e),(c)+hν+m+n,(d)−hν+m−n,(g)+hν+n,1+hν−(k)+n,(e)+m(j)+hν+n,1+hν−(h)ν′+n,(f)+m]xm​zhν+n​(−1)n​(H−K)m!​n!.\Sigma_{B}(z)=\sum_{\nu=1}^{B}\sum_{m=0}^{\infty}\sum_{n=0}^{\infty}\Gamma\Bigg{[}\begin{matrix}(b)_{\nu}^{\prime}-b_{\nu},&(h)-b_{\nu},&1+b_{\nu}-(b)_{\nu}^{\prime},&1+b_{\nu}-(h),&(f),\\ (d)-b_{\nu},&(k)-b_{\nu},&1+b_{\nu}-(d),&1+b_{\nu}-(k),&(e),\end{matrix}\\ \begin{matrix}(a)+b_{\nu}+2m+n,&(g)+b_{\nu}+m+n,&(e)+m,&1+b_{\nu}-(d)+n,&1+b_{\nu}-(k)+m+n\\ (c)+b_{\nu}+2m+n,&(j)+b_{\nu}+m+n,&(f)+m,&1+b_{\nu}-(b)_{\nu}^{\prime}+n,&1+b_{\nu}-(h)+m+n\end{matrix}\Bigg{]}\\ \times\frac{x^{m}z^{b_{\nu}+m+n}(-1)^{n(B+H-K-D)+m(H-K)}}{m!n!}\\ +\sum_{\nu=1}^{H}\sum_{m=0}^{\infty}\sum_{n=0}^{\infty}\Gamma\Bigg{[}\begin{matrix}(h)_{\nu}^{\prime}-h_{\nu},&1+h_{\nu}-(h)_{\nu}^{\prime},&(f),&(a)+h_{\nu}+m+n,&(b)-h_{\nu}+m-n,\\ (k)-h_{\nu},&1+h_{\nu}-(k),&(e),&(c)+h_{\nu}+m+n,&(d)-h_{\nu}+m-n,\end{matrix}\\ \begin{matrix}(g)+h_{\nu}+n,&1+h_{\nu}-(k)+n,&(e)+m\\ (j)+h_{\nu}+n,&1+h_{\nu}-(h)_{\nu}^{\prime}+n,&(f)+m\end{matrix}\Bigg{]}\frac{x^{m}z^{h_{\nu}+n}(-1)^{n(H-K)}}{m!n!}. (31)

Then after the integration over σ5\sigma_{5}, the result returns back to a strand of Γ\Gamma-function summations,

I=e3​ϵ​γ​(m42)−ν+3​D2Γ​(ν1)​Γ​(ν2)​Γ​(ν3)​Γ​(ν4)​Γ​(D/2)​∑m=0∞∑n=0∞1m!​n!​(m42m12)m​1(2​π​i)3​∭−i​∞i​∞dσ2​dσ3​dσ6×Γ​(−σ2)​Γ​(−σ3)​Γ​(−σ6)​Γ​(−σ2+σ6+ν3−D/2)​Γ​(σ2+σ3+D/2)×{(m42m12)σ2,3+ν1(m32m42)σ3(m22m42)−n+σ3−ν2+D2Γ[−σ2+σ6−ν1,2+ν−D,σ2−σ6+ν1,2,4−ν+D2+n,2​σ2+σ3−σ6−ν3+D+m+n,σ2,3+ν1+m,σ2−σ6+ν1,2−ν+D+12​σ2+σ3−σ6−ν3+D+n,σ3+ν2−D2+n,σ2−σ6+ν1,2−ν+D+1+m+n]+(m42m12)σ2,3+ν1​(m32m42)σ6​(m22m42)−n+σ2,3−σ6+ν1−ν+3​D2​Γ​[nm+n,ν4−D2+n]×Γ[σ2,3+ν4+m+n,σ2,3+ν1+m,−σ2+σ6−ν1,2+νD+1,σ2,3+ν4+n,−σ2+σ6−ν1,2+ν−D+1+n,σ2−σ6+ν1,2−ν+D−σ2,3+σ6−ν1+ν−3​D2+n]+(m42m12)ν−3​D2(m32m12)σ6(m22m12)nΓ[σ2,3−σ6+ν1,4−ν+D,2​σ2,3−σ6+ν1,4−ν+3​D2,σ2,3−σ6+ν1−ν+3​D2,−σ2,3+σ6−ν1,4+ν−D+1,−2​σ2,3+σ6−ν1,4+ν−3​D2+1,−σ2,3+σ6−ν1+ν−3​D2+1+m+n,−σ2,3+σ6−ν1+ν−3​D2+1,σ6+ν−3​D2+m+n,−σ2,3+σ6−ν1,4+ν−D+1+n,−2​σ2,3+σ6−ν1,4+ν−3​D2+1+n,σ2,3+ν4+m]+(m42m12)ν−3​D2​(m32m12)σ6​(m22m12)σ2,3−σ6+ν1,4−ν+D+n​Γ​[−ν4+D2,ν4−D2+1ν4−D2+1+m+n]×Γ[−σ2,3+σ6−ν1,4+ν−D,σ2,3+ν4+m,σ2,3−σ6+ν1,4−ν+D+1,σ2,3+D2,σ2,3+ν1,4−D2+m+n,σ2−σ6+ν1,2,4−ν+D2+n,−σ2,3−D2+1+n−σ2,3−D2+1,σ2,3−σ6+ν1,4−ν+D+1+n]+(m42m12)ν−3​D2​(m32m12)σ6​(m22m12)σ2,3−σ6+ν1−ν+3​D2+n​Γ​[ν4−D2,−ν4+D2+1,1+n1+m+n,−ν4+D2+1+n]Γ[−σ2,3+σ6−ν1+ν−3​D2,σ2,3−σ6+ν1−ν+3​D2+1,σ2,3+ν1+m+n,−σ2,3−ν4+1,σ2,3+ν4,−σ2,3−ν4+1+n,σ2,3+ν4+mσ2,3−σ6+ν1−ν+3​D2+1+n]},I=\frac{e^{3\epsilon\gamma}(m_{4}^{2})^{-\nu+\frac{3D}{2}}}{\Gamma(\nu_{1})\Gamma(\nu_{2})\Gamma(\nu_{3})\Gamma(\nu_{4})\Gamma(D/2)}\sum_{m=0}^{\infty}\sum_{n=0}^{\infty}\frac{1}{m!n!}\bigg{(}\frac{m_{4}^{2}}{m_{1}^{2}}\bigg{)}^{m}\frac{1}{(2\pi i)^{3}}\iiint_{-i\infty}^{i\infty}\mathrm{d}\sigma_{2}\mathrm{d}\sigma_{3}\mathrm{d}\sigma_{6}\\ \times\Gamma(-\sigma_{2})\Gamma(-\sigma_{3})\Gamma(-\sigma_{6})\Gamma(-\sigma_{2}+\sigma_{6}+\nu_{3}-D/2)\Gamma(\sigma_{2}+\sigma_{3}+D/2)\\ \times\Bigg{\{}\bigg{(}\frac{m_{4}^{2}}{m_{1}^{2}}\bigg{)}^{\sigma_{2,3}+\nu_{1}}\bigg{(}\frac{m_{3}^{2}}{m_{4}^{2}}\bigg{)}^{\sigma_{3}}\bigg{(}\frac{m_{2}^{2}}{m_{4}^{2}}\bigg{)}^{-n+\sigma_{3}-\nu_{2}+\frac{D}{2}}\Gamma\Bigg{[}\begin{matrix}-\sigma_{2}+\sigma_{6}-\nu_{1,2}+\nu-D,\\ \sigma_{2}-\sigma_{6}+\nu_{1,2,4}-\nu+\frac{D}{2}+n,\end{matrix}\\ \begin{matrix}2\sigma_{2}+\sigma_{3}-\sigma_{6}-\nu_{3}+D+m+n,&\sigma_{2,3}+\nu_{1}+m,&\sigma_{2}-\sigma_{6}+\nu_{1,2}-\nu+D+1\\ 2\sigma_{2}+\sigma_{3}-\sigma_{6}-\nu_{3}+D+n,&\sigma_{3}+\nu_{2}-\frac{D}{2}+n,&\sigma_{2}-\sigma_{6}+\nu_{1,2}-\nu+D+1+m+n\end{matrix}\Bigg{]}\\ +\bigg{(}\frac{m_{4}^{2}}{m_{1}^{2}}\bigg{)}^{\sigma_{2,3}+\nu_{1}}\bigg{(}\frac{m_{3}^{2}}{m_{4}^{2}}\bigg{)}^{\sigma_{6}}\bigg{(}\frac{m_{2}^{2}}{m_{4}^{2}}\bigg{)}^{-n+\sigma_{2,3}-\sigma_{6}+\nu_{1}-\nu+\frac{3D}{2}}\Gamma\Bigg{[}\begin{matrix}n\\ m+n,\ \nu_{4}-\frac{D}{2}+n\end{matrix}\Bigg{]}\\ \times\Gamma\Bigg{[}\begin{matrix}\sigma_{2,3}+\nu_{4}+m+n,&\sigma_{2,3}+\nu_{1}+m,&-\sigma_{2}+\sigma_{6}-\nu_{1,2}+\nu_{D}+1,\\ \sigma_{2,3}+\nu_{4}+n,&-\sigma_{2}+\sigma_{6}-\nu_{1,2}+\nu-D+1+n,&\end{matrix}\\ \begin{matrix}\sigma_{2}-\sigma_{6}+\nu_{1,2}-\nu+D\\ -\sigma_{2,3}+\sigma_{6}-\nu_{1}+\nu-\frac{3D}{2}+n\end{matrix}\Bigg{]}\\ +\bigg{(}\frac{m_{4}^{2}}{m_{1}^{2}}\bigg{)}^{\nu-\frac{3D}{2}}\bigg{(}\frac{m_{3}^{2}}{m_{1}^{2}}\bigg{)}^{\sigma_{6}}\bigg{(}\frac{m_{2}^{2}}{m_{1}^{2}}\bigg{)}^{n}\Gamma\Bigg{[}\begin{matrix}\sigma_{2,3}-\sigma_{6}+\nu_{1,4}-\nu+D,\\ 2\sigma_{2,3}-\sigma_{6}+\nu_{1,4}-\nu+\frac{3D}{2},\end{matrix}\\ \begin{matrix}\\ \sigma_{2,3}-\sigma_{6}+\nu_{1}-\nu+\frac{3D}{2},&-\sigma_{2,3}+\sigma_{6}-\nu_{1,4}+\nu-D+1,\\ -2\sigma_{2,3}+\sigma_{6}-\nu_{1,4}+\nu-\frac{3D}{2}+1,&-\sigma_{2,3}+\sigma_{6}-\nu_{1}+\nu-\frac{3D}{2}+1+m+n,\end{matrix}\\ \begin{matrix}-\sigma_{2,3}+\sigma_{6}-\nu_{1}+\nu-\frac{3D}{2}+1,&\sigma_{6}+\nu-\frac{3D}{2}+m+n,\\ -\sigma_{2,3}+\sigma_{6}-\nu_{1,4}+\nu-D+1+n,&\end{matrix}\\ \begin{matrix}-2\sigma_{2,3}+\sigma_{6}-\nu_{1,4}+\nu-\frac{3D}{2}+1+n,\ \sigma_{2,3}+\nu_{4}+m\\ \ \end{matrix}\Bigg{]}\\ +\bigg{(}\frac{m_{4}^{2}}{m_{1}^{2}}\bigg{)}^{\nu-\frac{3D}{2}}\bigg{(}\frac{m_{3}^{2}}{m_{1}^{2}}\bigg{)}^{\sigma_{6}}\bigg{(}\frac{m_{2}^{2}}{m_{1}^{2}}\bigg{)}^{\sigma_{2,3}-\sigma_{6}+\nu_{1,4}-\nu+D+n}\Gamma\Bigg{[}\begin{matrix}-\nu_{4}+\frac{D}{2},\ \nu_{4}-\frac{D}{2}+1\\ \nu_{4}-\frac{D}{2}+1+m+n\end{matrix}\Bigg{]}\\ \times\Gamma\Bigg{[}\begin{matrix}-\sigma_{2,3}+\sigma_{6}-\nu_{1,4}+\nu-D,\ \sigma_{2,3}+\nu_{4}+m,\ \sigma_{2,3}-\sigma_{6}+\nu_{1,4}-\nu+D+1,\\ \sigma_{2,3}+\frac{D}{2},\end{matrix}\\ \begin{matrix}\sigma_{2,3}+\nu_{1,4}-\frac{D}{2}+m+n,\ \sigma_{2}-\sigma_{6}+\nu_{1,2,4}-\nu+\frac{D}{2}+n,\ -\sigma_{2,3}-\frac{D}{2}+1+n\\ -\sigma_{2,3}-\frac{D}{2}+1,\ \sigma_{2,3}-\sigma_{6}+\nu_{1,4}-\nu+D+1+n\end{matrix}\Bigg{]}\\ +\bigg{(}\frac{m_{4}^{2}}{m_{1}^{2}}\bigg{)}^{\nu-\frac{3D}{2}}\bigg{(}\frac{m_{3}^{2}}{m_{1}^{2}}\bigg{)}^{\sigma_{6}}\bigg{(}\frac{m_{2}^{2}}{m_{1}^{2}}\bigg{)}^{\sigma_{2,3}-\sigma_{6}+\nu_{1}-\nu+\frac{3D}{2}+n}\Gamma\Bigg{[}\begin{matrix}\nu_{4}-\frac{D}{2},&-\nu_{4}+\frac{D}{2}+1,&1+n\\ 1+m+n,&-\nu_{4}+\frac{D}{2}+1+n&\end{matrix}\Bigg{]}\\ \Gamma\Bigg{[}\begin{matrix}-\sigma_{2,3}+\sigma_{6}-\nu_{1}+\nu-\frac{3D}{2},&\sigma_{2,3}-\sigma_{6}+\nu_{1}-\nu+\frac{3D}{2}+1,&\sigma_{2,3}+\nu_{1}+m+n,\\ -\sigma_{2,3}-\nu_{4}+1,&\sigma_{2,3}+\nu_{4},&\end{matrix}\\ \begin{matrix}-\sigma_{2,3}-\nu_{4}+1+n,\ \sigma_{2,3}+\nu_{4}+m\\ \sigma_{2,3}-\sigma_{6}+\nu_{1}-\nu+\frac{3D}{2}+1+n\end{matrix}\Bigg{]}\Bigg{\}}, (32)

where the Meijer theorem can be applied again. In the above equation, we have used the abbreviation σi,j,⋯=σi+σj+⋯\sigma_{i,j,\cdots}=\sigma_{i}+\sigma_{j}+\cdots.

Note that in some addends (i.e. the first and the third ones) the coefficients of σ2\sigma_{2} and σ3\sigma_{3} are not equal to 1 or -1. In this case, we can use Theorem 2.1 and imitate what we have done in eqs. (2.3)-(13). We have checked that the remaining three-fold integration can be done following such a strategy.

D6D_{6} diagrams with six internal legs

Another three-loop vacuum bubble worth discussion is the six-propagator D6D_{6} diagram  kniehl2017three ; Broadhurst:1998rz ; Freitas:2016zmy , as shown in figure 2.

k1−k2k_{1}-k_{2}k2−k3k_{2}-k_{3}k3−k1k_{3}-k_{1}k1k_{1}k3k_{3}k2k_{2}
Figure 2: The three-loop massive vacuum diagram D6D_{6}.

The definition of this integral is given by

D6≡e3​ϵ​γ​∫[d​k1]​[d​k2]​[d​k3]C1ν1​C2ν2​C3ν3​C4ν4​C5ν5​C6ν6,D_{6}\equiv e^{3\epsilon\gamma}\int\frac{[\mathrm{d}k_{1}][\mathrm{d}k_{2}][\mathrm{d}k_{3}]}{C_{1}^{\nu_{1}}C_{2}^{\nu_{2}}C_{3}^{\nu_{3}}C_{4}^{\nu_{4}}C_{5}^{\nu_{5}}C_{6}^{\nu_{6}}}, (33)

where

[d​ki]≡dD​kii​πD/2,i=1,2,3[\mathrm{d}k_{i}]\equiv\frac{\mathrm{d}^{D}k_{i}}{i\pi^{D/2}},\qquad i=1,2,3 (34)

and

Ci\displaystyle C_{i} ≡−ki2+mi2,i=1,2,3,\displaystyle\equiv-k_{i}^{2}+m_{i}^{2},\qquad i=1,2,3, (35)
C4\displaystyle C_{4} ≡−(k1−k2)2+m42,\displaystyle\equiv-(k_{1}-k_{2})^{2}+m_{4}^{2}, (36)
C5\displaystyle C_{5} ≡−(k2−k3)2+m52,\displaystyle\equiv-(k_{2}-k_{3})^{2}+m_{5}^{2}, (37)
C6\displaystyle C_{6} ≡−(k3−k1)2+m62.\displaystyle\equiv-(k_{3}-k_{1})^{2}+m_{6}^{2}. (38)

Its representations in Feynman parameters is

D6=e3​ϵ​γ​Γ​(ν−3​D2)Γ​(ν1)​Γ​(ν2)​Γ​(ν3)​Γ​(ν4)​Γ​(ν5)​Γ​(ν6)​∫01da1​da2​da3​da4​da5​da6δ​(1−∑i=16ai)​∏i=16aνi−1​𝒜−D2ℬν−3​D2,D_{6}=\frac{e^{3\epsilon\gamma}\Gamma(\nu-\frac{3D}{2})}{\Gamma(\nu_{1})\Gamma(\nu_{2})\Gamma(\nu_{3})\Gamma(\nu_{4})\Gamma(\nu_{5})\Gamma(\nu_{6})}\int_{0}^{1}\mathrm{d}a_{1}\mathrm{d}a_{2}\mathrm{d}a_{3}\mathrm{d}a_{4}\mathrm{d}a_{5}\mathrm{d}a_{6}\\ \delta\bigg{(}1-\sum_{i=1}^{6}a_{i}\bigg{)}\prod_{i=1}^{6}a^{\nu_{i}-1}\frac{\mathcal{A}^{-\frac{D}{2}}}{\mathcal{B}^{\nu-\frac{3D}{2}}}, (39)

where

𝒜=a1​a2​a3+a1​a3​a4+a2​a3​a4+a1​a2​a5+a1​a3​a5+a1​a4​a5+a2​a4​a5+a1​a3​a6+a2​a3​a6+a1​a4​a6+a2​a4​a6+a3​a4​a6+a1​a5​a6+a2​a5​a6+a3​a5​a6ℬ=a1​m12+a2​m22+a3​m32+a4​m42+a5​m52+a6​m62.\displaystyle\begin{split}\mathcal{A}&=a_{1}a_{2}a_{3}+a_{1}a_{3}a_{4}+a_{2}a_{3}a_{4}+a_{1}a_{2}a_{5}+a_{1}a_{3}a_{5}+a_{1}a_{4}a_{5}+a_{2}a_{4}a_{5}+a_{1}a_{3}a_{6}\\ &+a_{2}a_{3}a_{6}+a_{1}a_{4}a_{6}+a_{2}a_{4}a_{6}+a_{3}a_{4}a_{6}+a_{1}a_{5}a_{6}+a_{2}a_{5}a_{6}+a_{3}a_{5}a_{6}\\ \mathcal{B}&=a_{1}m_{1}^{2}+a_{2}m_{2}^{2}+a_{3}m_{3}^{2}+a_{4}m_{4}^{2}+a_{5}m_{5}^{2}+a_{6}m_{6}^{2}.\end{split} (40)

After introducing the Mellin-Barnes transformation to both the numerator and denominator, we have

D6=e3​ϵ​γ​Γ​(ν−3​D2)Γ​(ν1)​Γ​(ν2)​Γ​(ν3)​Γ​(ν4)​Γ​(ν5)​Γ​(ν6)​∫01da1​da2​da3​da4​da5​da6​δ​(1−∑i=16ai)×1(2​π​i)20​∫−i​∞i​∞dσ1​⋯​dσ20​∏i=120Γ​(−σi)​Γ​(∑i=115σi+D2)​Γ​(∑i=1620σi+ν−3​D2)×(m12m62)σ16​(m22m62)σ17​(m32m62)σ18​(m42m62)σ19​(m52m62)σ20​(m62)−ν+3​D2×a1σ1,2,4,5,6,9,11,14,16+ν1−1​a2σ1,3,4,7,9,10,12,15,17+ν2−1​a3−σ4,6,7,9,11,12,14,15+σ18−D2+ν3−1×a4σ2,3,6,7,8,11,12,13,19+ν4−1​a5−σ1,2,3,9,10,11,12,13+σ20−D2+ν5−1​a6−σ1,2,3,4,5,6,7,8,16,17,18,19,20−ν+D+ν6−1.D_{6}=\frac{e^{3\epsilon\gamma}\Gamma(\nu-\frac{3D}{2})}{\Gamma(\nu_{1})\Gamma(\nu_{2})\Gamma(\nu_{3})\Gamma(\nu_{4})\Gamma(\nu_{5})\Gamma(\nu_{6})}\int_{0}^{1}\mathrm{d}a_{1}\mathrm{d}a_{2}\mathrm{d}a_{3}\mathrm{d}a_{4}\mathrm{d}a_{5}\mathrm{d}a_{6}\delta\bigg{(}1-\sum_{i=1}^{6}a_{i}\bigg{)}\\ \times\frac{1}{(2\pi i)^{20}}\int_{-i\infty}^{i\infty}\mathrm{d}\sigma_{1}\cdots\mathrm{d}\sigma_{20}\prod_{i=1}^{20}\Gamma(-\sigma_{i})\Gamma\bigg{(}\sum_{i=1}^{15}\sigma_{i}+\frac{D}{2}\bigg{)}\Gamma\bigg{(}\sum_{i=16}^{20}\sigma_{i}+\nu-\frac{3D}{2}\bigg{)}\\ \times\bigg{(}\frac{m_{1}^{2}}{m_{6}^{2}}\bigg{)}^{\sigma_{16}}\bigg{(}\frac{m_{2}^{2}}{m_{6}^{2}}\bigg{)}^{\sigma_{17}}\bigg{(}\frac{m_{3}^{2}}{m_{6}^{2}}\bigg{)}^{\sigma_{18}}\bigg{(}\frac{m_{4}^{2}}{m_{6}^{2}}\bigg{)}^{\sigma_{19}}\bigg{(}\frac{m_{5}^{2}}{m_{6}^{2}}\bigg{)}^{\sigma_{20}}(m_{6}^{2})^{-\nu+\frac{3D}{2}}\\ \times a_{1}^{\sigma_{1,2,4,5,6,9,11,14,16}+\nu_{1}-1}a_{2}^{\sigma_{1,3,4,7,9,10,12,15,17}+\nu_{2}-1}a_{3}^{-\sigma_{4,6,7,9,11,12,14,15}+\sigma_{18}-\frac{D}{2}+\nu_{3}-1}\\ \times a_{4}^{\sigma_{2,3,6,7,8,11,12,13,19}+\nu_{4}-1}a_{5}^{-\sigma_{1,2,3,9,10,11,12,13}+\sigma_{20}-\frac{D}{2}+\nu_{5}-1}a_{6}^{-\sigma_{1,2,3,4,5,6,7,8,16,17,18,19,20}-\nu+D+\nu_{6}-1}. (41)

Integrating over Feynman parameters and setting mi=1m_{i}=1, we obtain the Mellin-Barnes integral,

D6=e3​ϵ​γ​Γ​(ν−3​D2)Γ​(ν1)​Γ​(ν2)​Γ​(ν3)​Γ​(ν4)​Γ​(ν5)​Γ​(ν6)​1(2​π​i)20​∫−i​∞i​∞dσ1​⋯​dσ20​∏i=120Γ​(−σi)×Γ​(∑i=115σi+D2)​Γ​(∑i=1620σi+ν−3​D2)​Γ​(σ1,2,4,5,6,9,11,14,16+ν1)×Γ​(σ1,3,4,7,9,10,12,15,17+ν2)​Γ​(−σ4,6,7,9,11,12,14,15+σ18−D2+ν3)×Γ​(σ2,3,6,7,8,11,12,13,19+ν4)​Γ​(−σ1,2,3,9,10,11,12,13+σ20−D2+ν5)×Γ​(−σ1,2,3,4,5,6,7,8,16,17,18,19,20−ν+D+ν6).D_{6}=\frac{e^{3\epsilon\gamma}\Gamma(\nu-\frac{3D}{2})}{\Gamma(\nu_{1})\Gamma(\nu_{2})\Gamma(\nu_{3})\Gamma(\nu_{4})\Gamma(\nu_{5})\Gamma(\nu_{6})}\frac{1}{(2\pi i)^{20}}\int_{-i\infty}^{i\infty}\mathrm{d}\sigma_{1}\cdots\mathrm{d}\sigma_{20}\prod_{i=1}^{20}\Gamma(-\sigma_{i})\\ \times\Gamma\bigg{(}\sum_{i=1}^{15}\sigma_{i}+\frac{D}{2}\bigg{)}\Gamma\bigg{(}\sum_{i=16}^{20}\sigma_{i}+\nu-\frac{3D}{2}\bigg{)}\Gamma(\sigma_{1,2,4,5,6,9,11,14,16}+\nu_{1})\\ \times\Gamma(\sigma_{1,3,4,7,9,10,12,15,17}+\nu_{2})\Gamma(-\sigma_{4,6,7,9,11,12,14,15}+\sigma_{18}-\frac{D}{2}+\nu_{3})\\ \times\Gamma(\sigma_{2,3,6,7,8,11,12,13,19}+\nu_{4})\Gamma(-\sigma_{1,2,3,9,10,11,12,13}+\sigma_{20}-\frac{D}{2}+\nu_{5})\\ \times\Gamma(-\sigma_{1,2,3,4,5,6,7,8,16,17,18,19,20}-\nu+D+\nu_{6}). (42)

The form of the integrand becomes obviously the product of Γ\Gamma functions and exponential functions again, though there are 20 Mellin-Barnes parameters. Hence this calculation can be done by applying the Meijer theorem and its corollary in turns. For instance, if we integrate over σ1\sigma_{1}, we obtain

D6=e3​ϵ​γ​Γ​(ν−3​D2)Γ​(ν1)​Γ​(ν2)​Γ​(ν3)​Γ​(ν4)​Γ​(ν5)​Γ​(ν6)​1(2​π​i)19​∫−i​∞i​∞dσ2​⋯​dσ20​∏i=220Γ​(−σi)×Γ​(∑i=1620σi+ν−3​D2)​Γ​(−σ4,6,7,9,11,12,14,15+σ18−D2+ν3)​Γ​(σ2,3,6,7,8,11,12,13,19+ν4)×{Γ(∑i=215σi+D2)Γ(σ2,4,5,6,9,11,14,16+ν1)Γ(σ3,4,7,9,10,12,15,17+ν2)×Γ​(−σ2,3,9,10,11,12,13+σ20−D2+ν5)​Γ​(−σ2,3,4,5,6,7,8,16,17,18,19,20−ν+D+ν6)×F23​[𝒜1,𝒜2,𝒜3;𝒜4,𝒜5;−1]+Γ​(σ4,5,6,7,8,14,15,20+ν5)​Γ​(−σ3,10,12,13+σ4,5,6,14,16,20−D2+ν1+ν5)×Γ​(−σ2,11,13+σ4,7,15,17,20−D2+ν2+ν5)​Γ​(σ2,3,9,10,11,12,13−σ20+D2−ν5)×Γ​(σ4,5,6,7,8,16,17,19−σ9,10,11,12,13+2​σ20+ν−3​D2−ν6+ν5)×F23​[ℬ1,ℬ2,ℬ3;ℬ4,ℬ5;−1]+Γ​(σ9,10,11,12,13,14,15−σ16,17,18,19,20−ν+3​D2+ν6)​Γ​(−σ3,7,8,17,18,19,20+σ9,11,14−ν+D+ν1+ν6)×Γ​(−σ2,5,6,8,16,18,19,20+σ9,10,12,15−ν+D+ν2+ν6)​Γ​(σ2,3,4,5,6,7,8,,16,17,18,19,20+ν−D−ν6)×Γ​(σ4,5,6,7,8,16,17,18,19−σ9,10,11,12,13+2​σ20+ν−3​D2+ν5−ν6)×F23[𝒞1,𝒞2,𝒞3;𝒞4,𝒞5;−1]},D_{6}=\frac{e^{3\epsilon\gamma}\Gamma(\nu-\frac{3D}{2})}{\Gamma(\nu_{1})\Gamma(\nu_{2})\Gamma(\nu_{3})\Gamma(\nu_{4})\Gamma(\nu_{5})\Gamma(\nu_{6})}\frac{1}{(2\pi i)^{19}}\int_{-i\infty}^{i\infty}\mathrm{d}\sigma_{2}\cdots\mathrm{d}\sigma_{20}\prod_{i=2}^{20}\Gamma(-\sigma_{i})\\ \times\Gamma\bigg{(}\sum_{i=16}^{20}\sigma_{i}+\nu-\frac{3D}{2}\bigg{)}\Gamma(-\sigma_{4,6,7,9,11,12,14,15}+\sigma_{18}-\frac{D}{2}+\nu_{3})\Gamma(\sigma_{2,3,6,7,8,11,12,13,19}+\nu_{4})\\ \times\Bigg{\{}\Gamma\bigg{(}\sum_{i=2}^{15}\sigma_{i}+\frac{D}{2}\bigg{)}\Gamma(\sigma_{2,4,5,6,9,11,14,16}+\nu_{1})\Gamma(\sigma_{3,4,7,9,10,12,15,17}+\nu_{2})\\ \times\Gamma(-\sigma_{2,3,9,10,11,12,13}+\sigma_{20}-\frac{D}{2}+\nu_{5})\Gamma(-\sigma_{2,3,4,5,6,7,8,16,17,18,19,20}-\nu+D+\nu_{6})\\ \times{}_{3}F_{2}\bigg{[}\begin{matrix}\mathcal{A}_{1},&\mathcal{A}_{2},&\mathcal{A}_{3};\\ \mathcal{A}_{4},&\mathcal{A}_{5}&;\end{matrix}\quad-1\bigg{]}\\ +\Gamma(\sigma_{4,5,6,7,8,14,15,20}+\nu_{5})\Gamma(-\sigma_{3,10,12,13}+\sigma_{4,5,6,14,16,20}-\frac{D}{2}+\nu_{1}+\nu_{5})\\ \times\Gamma(-\sigma_{2,11,13}+\sigma_{4,7,15,17,20}-\frac{D}{2}+\nu_{2}+\nu_{5})\Gamma(\sigma_{2,3,9,10,11,12,13}-\sigma_{20}+\frac{D}{2}-\nu_{5})\\ \times\Gamma(\sigma_{4,5,6,7,8,16,17,19}-\sigma_{9,10,11,12,13}+2\sigma_{20}+\nu-\frac{3D}{2}-\nu_{6}+\nu_{5})\\ \times{}_{3}F_{2}\bigg{[}\begin{matrix}\mathcal{B}_{1},&\mathcal{B}_{2},&\mathcal{B}_{3};\\ \mathcal{B}_{4},&\mathcal{B}_{5}&;\end{matrix}\quad-1\bigg{]}\\ +\Gamma(\sigma_{9,10,11,12,13,14,15}-\sigma_{16,17,18,19,20}-\nu+\frac{3D}{2}+\nu_{6})\Gamma(-\sigma_{3,7,8,17,18,19,20}+\sigma_{9,11,14}-\nu+D+\nu_{1}+\nu_{6})\\ \times\Gamma(-\sigma_{2,5,6,8,16,18,19,20}+\sigma_{9,10,12,15}-\nu+D+\nu_{2}+\nu_{6})\Gamma(\sigma_{2,3,4,5,6,7,8,,16,17,18,19,20}+\nu-D-\nu_{6})\\ \times\Gamma(\sigma_{4,5,6,7,8,16,17,18,19}-\sigma_{9,10,11,12,13}+2\sigma_{20}+\nu-\frac{3D}{2}+\nu_{5}-\nu_{6})\\ \times{}_{3}F_{2}\bigg{[}\begin{matrix}\mathcal{C}_{1},&\mathcal{C}_{2},&\mathcal{C}_{3};\\ \mathcal{C}_{4},&\mathcal{C}_{5}&;\end{matrix}\quad-1\bigg{]}\Bigg{\}}, (43)

where the explicit expressions of 𝒜i,ℬi\mathcal{A}_{i},\mathcal{B}_{i}, and 𝒞i\mathcal{C}_{i} are collected in appendix B.

The following computation becomes too tedious to be shown explicitly here, though we do not expect any difficulties in principle.

ppk1k_{1}k2k_{2}p−k1−k2p-k_{1}-k_{2}pp
Figure 3: The two-loop sunset diagram.

Diagrams with external legs

Now, let us extend our exploration to the cases of diagrams with external legs. Take the well-known sunset diagram in figure 3 as an example. The Feynman integral is defined as

S≡e2​ϵ​γ​∬[d​k1]​[d​k2]C1ν1​C2ν2​C3ν3,S\equiv e^{2\epsilon\gamma}\iint\frac{[\mathrm{d}k_{1}][\mathrm{d}k_{2}]}{C_{1}^{\nu_{1}}C_{2}^{\nu_{2}}C_{3}^{\nu_{3}}}, (44)

where

[d​ki]\displaystyle[\mathrm{d}k_{i}] ≡dD​kii​πD/2,i=1,2,\displaystyle\equiv\frac{\mathrm{d}^{D}k_{i}}{i\pi^{D/2}},\qquad i=1,2, (45)
Ci\displaystyle C_{i} ≡−ki2+mi2,i=1,2,\displaystyle\equiv-k_{i}^{2}+m_{i}^{2},\qquad i=1,2, (46)
C3\displaystyle C_{3} ≡−(−k1−k2+p)2+m32.\displaystyle\equiv-(-k_{1}-k_{2}+p)^{2}+m_{3}^{2}. (47)

Adopting the Feynman parameterization, we have the following form,

S=e2​ϵ​γ​Γ​(ν−D)Γ​(ν1)​Γ​(ν2)​Γ​(ν3)​∫01da1​da2​da3​δ​(1−a1−a2−a3)​a1ν1−1​a2ν2−1​a3ν3−1×(a1​a2+a1​a3+a2​a3)ν−3​D2[(a1​a2+a1​a3+a2​a3)​(a1​m12+a2​m22+a3​m32−a3​p2)+a32​p2​(a1+a2)]ν−D.S=\frac{e^{2\epsilon\gamma}\Gamma(\nu-D)}{\Gamma(\nu_{1})\Gamma(\nu_{2})\Gamma(\nu_{3})}\int_{0}^{1}\mathrm{d}a_{1}\mathrm{d}a_{2}\mathrm{d}a_{3}\delta(1-a_{1}-a_{2}-a_{3})a_{1}^{\nu_{1}-1}a_{2}^{\nu_{2}-1}a_{3}^{\nu_{3}-1}\\ \times\frac{(a_{1}a_{2}+a_{1}a_{3}+a_{2}a_{3})^{\nu-\frac{3D}{2}}}{[(a_{1}a_{2}+a_{1}a_{3}+a_{2}a_{3})(a_{1}m_{1}^{2}+a_{2}m_{2}^{2}+a_{3}m_{3}^{2}-a_{3}p^{2})+a_{3}^{2}p^{2}(a_{1}+a_{2})]^{\nu-D}}. (48)

Its corresponding Mellin-Barnes transformation is given by

S=e2​ϵ​γΓ​(ν1)​Γ​(ν2)​Γ​(ν3)​Γ​(−ν+3​D2)​1(2​π​i)8​∫−i​∞i​∞dσ1​dσ2​dσ3​dσ4​dσ5​dσ6​dσ7​dσ8×Γ​(−σ1)​Γ​(−σ2)​Γ​(−σ3)​Γ​(−σ4)​Γ​(−σ5)​Γ​(−σ6)​Γ​(−σ7)​Γ​(−σ8)​Γ​(σ1,2−ν+3​D2)×Γ​(σ3,4,5,6,7,8+ν−D)​(m12)σ3,5​(m22)σ4,7​(m32)σ6,8​(m12+m22+m32−p2)−σ3,4,5,6,7,8−ν+D×Γ​(σ1,2,3,5−σ7,8−ν+D+ν1)​Γ​(−σ2,5,6+σ4,7+ν2−D2)​Γ​(−σ1,3,4+σ6,8+ν3−D2).S=\frac{e^{2\epsilon\gamma}}{\Gamma(\nu_{1})\Gamma(\nu_{2})\Gamma(\nu_{3})\Gamma(-\nu+\frac{3D}{2})}\frac{1}{(2\pi i)^{8}}\int_{-i\infty}^{i\infty}\mathrm{d}\sigma_{1}\mathrm{d}\sigma_{2}\mathrm{d}\sigma_{3}\mathrm{d}\sigma_{4}\mathrm{d}\sigma_{5}\mathrm{d}\sigma_{6}\mathrm{d}\sigma_{7}\mathrm{d}\sigma_{8}\\ \times\Gamma(-\sigma_{1})\Gamma(-\sigma_{2})\Gamma(-\sigma_{3})\Gamma(-\sigma_{4})\Gamma(-\sigma_{5})\Gamma(-\sigma_{6})\Gamma(-\sigma_{7})\Gamma(-\sigma_{8})\Gamma(\sigma_{1,2}-\nu+\frac{3D}{2})\\ \times\Gamma(\sigma_{3,4,5,6,7,8}+\nu-D)(m_{1}^{2})^{\sigma_{3,5}}(m_{2}^{2})^{\sigma_{4,7}}(m_{3}^{2})^{\sigma_{6,8}}(m_{1}^{2}+m_{2}^{2}+m_{3}^{2}-p^{2})^{-\sigma_{3,4,5,6,7,8}-\nu+D}\\ \times\Gamma(\sigma_{1,2,3,5}-\sigma_{7,8}-\nu+D+\nu_{1})\Gamma(-\sigma_{2,5,6}+\sigma_{4,7}+\nu_{2}-\frac{D}{2})\Gamma(-\sigma_{1,3,4}+\sigma_{6,8}+\nu_{3}-\frac{D}{2}). (49)

Again, we see that the integral is in the form ready for applying the Meijer theorem and its corollary. For example, let us perform the integration over σ1\sigma_{1} to demonstrate the feasibility of our algorithm, and the result reads

S=e2​ϵ​γΓ​(ν1)​Γ​(ν2)​Γ​(ν3)​Γ​(−ν+3​D2)​1(2​π​i)7​∫−i​∞i​∞dσ2​dσ3​dσ4​dσ5​dσ6​dσ7​dσ8×Γ​(−σ2)​Γ​(−σ3)​Γ​(−σ4)​Γ​(−σ5)​Γ​(−σ6)​Γ​(−σ7)​Γ​(−σ8)×Γ​(σ3,4,5,6,7,8+ν−D)​(m12)σ3,5​(m22)σ4,7​(m32)σ6,8(m12+m22+m32−p2)−σ3,4,5,6,7,8−ν+D​Γ​(−σ2,5,6+σ4,7+ν2−D2)​Γ​(σ2−ν+3​D2)×Γ​(σ2,6,8−σ3,4−ν+ν3+D)​Γ​(σ2,3,5−σ7,8−ν+D+ν1)​Γ​(σ2,5,6−σ4,7−ν2+D2)Γ​(2​σ2−σ4,7+σ5,6−ν−ν2+2​D).S=\frac{e^{2\epsilon\gamma}}{\Gamma(\nu_{1})\Gamma(\nu_{2})\Gamma(\nu_{3})\Gamma(-\nu+\frac{3D}{2})}\frac{1}{(2\pi i)^{7}}\int_{-i\infty}^{i\infty}\mathrm{d}\sigma_{2}\mathrm{d}\sigma_{3}\mathrm{d}\sigma_{4}\mathrm{d}\sigma_{5}\mathrm{d}\sigma_{6}\mathrm{d}\sigma_{7}\mathrm{d}\sigma_{8}\\ \times\Gamma(-\sigma_{2})\Gamma(-\sigma_{3})\Gamma(-\sigma_{4})\Gamma(-\sigma_{5})\Gamma(-\sigma_{6})\Gamma(-\sigma_{7})\Gamma(-\sigma_{8})\\ \times\Gamma(\sigma_{3,4,5,6,7,8}+\nu-D)(m_{1}^{2})^{\sigma_{3,5}}(m_{2}^{2})^{\sigma_{4,7}}(m_{3}^{2})^{\sigma_{6,8}}\\ (m_{1}^{2}+m_{2}^{2}+m_{3}^{2}-p^{2})^{-\sigma_{3,4,5,6,7,8}-\nu+D}\Gamma(-\sigma_{2,5,6}+\sigma_{4,7}+\nu_{2}-\frac{D}{2})\Gamma(\sigma_{2}-\nu+\frac{3D}{2})\\ \times\Gamma(\sigma_{2,6,8}-\sigma_{3,4}-\nu+\nu_{3}+D)\frac{\Gamma(\sigma_{2,3,5}-\sigma_{7,8}-\nu+D+\nu_{1})\Gamma(\sigma_{2,5,6}-\sigma_{4,7}-\nu_{2}+\frac{D}{2})}{\Gamma(2\sigma_{2}-\sigma_{4,7}+\sigma_{5,6}-\nu-\nu_{2}+2D)}. (50)

We do not bother to show the subsequent steps.

The calculation algorithm

Summarizing the previous instances, we are ready to lay down a calculation algorithm for Feynman integral calculations. Starting with a given Feynman integral, we transform all the integral variables into Feynman parameters, and then to Mellin-Barnes forms. After these steps, the integrand would appear to be in the form of multiplication of Γ\Gamma functions and exponential functions, and it is natural to apply the Meijer theorem and express the result in products of Γ\Gamma functions, exponential functions, and hypergeometric functions. The integration of these products can be dealt with by Corollary 2.1 and the result becomes integrals of Γ\Gamma functions and exponential functions again. The procedure will continue until the computation is eventually done. Note that before applying the theorems, it is essential to normalize the factors in front of the integration variables to ±1\pm 1 by using Theorem 2.1. A workflow is displayed in figure 4.

Feynman IntegralsMellin-Barnes RepresentationWhether the calculation is done Integration over Γ\Gamma functions and exponential functions Integration over Γ\Gamma functions, exponential functions and hypergeometric functions Analytical ResultsAsymptotic Series of ϵ\epsilonNumerical ResultsNOCorollary 2.1 & NormalizationYESMeijer Theorem & NormalizationAsymptotic Expansion in ϵ\epsilon
Figure 4: The workflow of the calculation of Feynman integrals with the Meijer theorem and its corollary.

4 Conclusion and outlook

We have proposed an analytical calculation method for the Feynman integrals in arbitrary dimensions by applying the Meijer theorem and its corollary, for which we have presented an explicit proof. The feasibility has been demonstrated by the calculation of the three-loop vacuum banana diagram and discussed for several other multi-loop Feynman diagrams. Even though our algorithm is general, there are still some obstacles arising in practical calculations. Here we list two of them.

  • •

    High computation complexity in multi-fold integration. It has been evident that in our discussion in section 3 we completed only the first few steps to show the feasibility of our algorithm. The expression starts to explode in the following steps. An automatic computer program encoding the algorithm is needed.

  • •

    Necessity of normalization. The Meijer theorem and its corollary only apply to the integral in which the argument of the Γ\Gamma function is normalized. In section 2, we have illustrated the normalization procedure using the theorem 2.1. However, one may encounter more complicated cases, such as Γ​(σ3+σ2)​Γ​(3​σ3+2​σ2)\Gamma(\sigma_{3}+\sigma_{2})\Gamma(3\sigma_{3}+2\sigma_{2}) in the integrand. The normalization of the σ2\sigma_{2} integration makes the σ3\sigma_{3} integration abnormal. Our method does not work in such cases.

Acknowledgements.
This work was supported in part by the National Science Foundation of China (grant Nos. 12005117, 12321005, 12375076) and the Taishan Scholar Foundation of Shandong province (tsqn201909011). The Feynman diagrams in this paper were drawn using the TikZ-Feynman package Ellis:2016jkw .

Appendix A The proof of Corollary 2.1

We begin with a careful consideration of the following complex integral

ℐ(z)=12​π​i∫−i​∞i​∞Γ[(a)+s,(b)−s,(g)+s,(h)−s(c)+s,(d)−s,(j)+s,(k)−s]×FC+D+FA+B+E​[(a)+s,(b)−s,(e);(c)+s,(d)−s,(f);x]zs​d​s.\mathcal{I}(z)=\frac{1}{2\pi i}\int_{-i\infty}^{i\infty}\Gamma\Bigg{[}\begin{matrix}(a)+s,&(b)-s,&(g)+s,(h)-s\\ (c)+s,&(d)-s,&(j)+s,(k)-s\end{matrix}\Bigg{]}\times\\ {}_{A+B+E}F_{C+D+F}\Bigg{[}\begin{matrix}(a)+s,&(b)-s,&(e);\\ (c)+s,&(d)-s,&(f);\end{matrix}\quad x\Bigg{]}\quad z^{s}\mathrm{d}s. (51)

The generalized hypergeometric functions can be expanded as

FC+D+FA+B+E​[(a)+s,(b)−s,(e);(c)+s,(d)−s,(f);x]≡∑m=0∞((a)+s)m​((b)−s)m​((e))m((c)+s)m​((d)−s)m​((f))m​xmm!,{}_{A+B+E}F_{C+D+F}\Bigg{[}\begin{matrix}(a)+s,&(b)-s,&(e);\\ (c)+s,&(d)-s,&(f);\end{matrix}\quad x\Bigg{]}\equiv\sum_{m=0}^{\infty}\frac{((a)+s)_{m}((b)-s)_{m}((e))_{m}}{((c)+s)_{m}((d)-s)_{m}((f))_{m}}\frac{x^{m}}{m!}, (52)

or in terms of Γ\Gamma functions,

FC+D+FA+B+E​[(a)+s,(b)−s,(e);(c)+s,(d)−s,(f);x]=∑m=0∞Γ​[(a)+s+m,(b)−s+m,(e)+m,(c)+s,(d)−s,(f)(a)+s,(b)−s,(e),(c)+s+m,(d)−s+m,(f)+m]​xmm!.{}_{A+B+E}F_{C+D+F}\Bigg{[}\begin{matrix}(a)+s,&(b)-s,&(e);\\ (c)+s,&(d)-s,&(f);\end{matrix}\quad x\Bigg{]}\\ =\sum_{m=0}^{\infty}\Gamma\Bigg{[}\begin{matrix}(a)+s+m,&(b)-s+m,&(e)+m,&(c)+s,&(d)-s,&(f)\\ (a)+s,&(b)-s,&(e),&(c)+s+m,&(d)-s+m,&(f)+m\end{matrix}\Bigg{]}\frac{x^{m}}{m!}. (53)

Then the previous integral (51) becomes

ℐ​(z)=∑m=0∞Γ​[(e)+m,(f)(e),(f)+m]​xmm!​12​π​i​∫−i​∞i​∞ds​zs×Γ​[(a)+s+m,(b)−s+m,(g)+s,(h)−s(c)+s+m,(d)−s+m,(j)+s,(k)−s].\mathcal{I}(z)=\sum_{m=0}^{\infty}\Gamma\Bigg{[}\begin{matrix}(e)+m,&(f)\\ (e),&(f)+m\end{matrix}\Bigg{]}\frac{x^{m}}{m!}\frac{1}{2\pi i}\int_{-i\infty}^{i\infty}\mathrm{d}sz^{s}\\ \times\Gamma\Bigg{[}\begin{matrix}(a)+s+m,&(b)-s+m,&(g)+s,&(h)-s\\ (c)+s+m,&(d)-s+m,&(j)+s,&(k)-s\end{matrix}\Bigg{]}\,. (54)

According to the Meijer theorem, the result would be

  1. (i)

    ℐ​(z)=ΣA​(z)\mathcal{I}(z)=\Sigma_{A}(z) when A+G+D+K>B+H+C+J\quad A+G+D+K>B+H+C+J\qquad

    or A+G+D+K=B+H+C+J\quad A+G+D+K=B+H+C+J and |z|>1|z|>1,

  2. (ii)

    ℐ​(z)=ΣB​(z)\mathcal{I}(z)=\Sigma_{B}(z) when A+G+D+K<B+H+C+J\quad A+G+D+K<B+H+C+J\qquad

    or A+G+D+K=B+H+C+J\quad A+G+D+K=B+H+C+J and |z|<1|z|<1,

provided that 12​π​|A+G+B+H−C−D−J−K|>|arg⁡z|\frac{1}{2}\pi|A+G+B+H-C-D-J-K|>|\arg z|. In the special case of z=1z=1 and Re​(Σ​c+Σ​d+Σ​j+Σ​k−Σ​a−Σ​b−Σ​g−Σ​h)>0\textrm{Re}(\Sigma c+\Sigma d+\Sigma j+\Sigma k-\Sigma a-\Sigma b-\Sigma g-\Sigma h)>0,

  1. (iii)

    ℐ​(1)=ΣA​(1)=ΣB​(1)\mathcal{I}(1)=\Sigma_{A}(1)=\Sigma_{B}(1) when A+G−C−J=B+H−D−K≥0\quad A+G-C-J=B+H-D-K\geq 0.

The explicit expressions of ΣA\Sigma_{A} and ΣB\Sigma_{B} are given by

ΣA​(z)=∑m=0∞Γ​[(e)+m,(f)(e),(f)+m]​xmm!×{∑μ=1Az−aμ−mΓ[(a)μ′−aμ,(g)−aμ−m,(b)+aμ+2​m,(h)+aμ+m(c)−aμ,(j)−aμ−m,(d)+aμ+2​m,(k)+aμ+m]×FA+G+D+K−1B+H+C+J[(b)+aμ+2​m,(h)+aμ+m,1+aμ−(c),1+aμ+m−(j);1+aμ−(a)μ′,1+aμ+m−(g),(d)+aμ+2​m,(k)+aμ+m;(−1)A+G+C+Jz−1]+∑μ=1Gz−gμΓ[(a)−gμ+m,(g)μ′−gμ,(b)+gμ+m(h)+gμ(c)−gμ+m,(j)−gμ,(d)+gμ+m,(k)+gμ]×FA+G+D+K−1B+H+C+J[(b)+gμ+m,(h)+gμ,1+gμ−(c)−m,1+gμ−(j);1+gμ−(a)−m,1+gμ−(g)μ′,(d)+gμ+m,(k)+gμ;(−1)A+G+C+Jz−1]}\Sigma_{A}(z)=\sum_{m=0}^{\infty}\Gamma\Bigg{[}\begin{matrix}(e)+m,&(f)\\ (e),&(f)+m\end{matrix}\Bigg{]}\frac{x^{m}}{m!}\\ \times\Bigg{\{}\sum_{\mu=1}^{A}z^{-a_{\mu}-m}\Gamma\Bigg{[}\begin{matrix}(a)_{\mu}^{\prime}-a_{\mu},&(g)-a_{\mu}-m,&(b)+a_{\mu}+2m,&(h)+a_{\mu}+m\\ (c)-a_{\mu},&(j)-a_{\mu}-m,&(d)+a_{\mu}+2m,&(k)+a_{\mu}+m\end{matrix}\Bigg{]}\times\\ {}_{B+H+C+J}F_{A+G+D+K-1}\Bigg{[}\begin{matrix}(b)+a_{\mu}+2m,&(h)+a_{\mu}+m,&1+a_{\mu}-(c),&1+a_{\mu}+m-(j);\\ 1+a_{\mu}-(a)_{\mu}^{\prime},&1+a_{\mu}+m-(g),&(d)+a_{\mu}+2m,&(k)+a_{\mu}+m;\end{matrix}\quad\\ (-1)^{A+G+C+J}z^{-1}\Bigg{]}\\ +\sum_{\mu=1}^{G}z^{-g_{\mu}}\Gamma\Bigg{[}\begin{matrix}(a)-g_{\mu}+m,&(g)_{\mu}^{\prime}-g_{\mu},&(b)+g_{\mu}+m&(h)+g_{\mu}\\ (c)-g_{\mu}+m,&(j)-g_{\mu},&(d)+g_{\mu}+m,&(k)+g_{\mu}\end{matrix}\Bigg{]}\times\\ {}_{B+H+C+J}F_{A+G+D+K-1}\Bigg{[}\begin{matrix}(b)+g_{\mu}+m,&(h)+g_{\mu},&1+g_{\mu}-(c)-m,&1+g_{\mu}-(j);\\ 1+g_{\mu}-(a)-m,&1+g_{\mu}-(g)_{\mu}^{\prime},&(d)+g_{\mu}+m,&(k)+g_{\mu};\end{matrix}\\ \quad(-1)^{A+G+C+J}z^{-1}\Bigg{]}\Bigg{\}} (55)

and

ΣB​(z)=∑m=0∞Γ​[(e)+m,(f)(e),(f)+m]​xmm!×{∑ν=1Bzbν+mΓ[(a)+bν+2​m,(g)+bν+m,(b)ν′−bν,(h)−bν−m(c)+bν+2​m,(j)+bν+m,(d)−bν,(k)−bν−m]×FB+H+C+J−1A+G+D+K[(a)+bν+2​m,(g)+bν+m,1+bν−(d),1+bν+m−(k);1+bν−(b)ν′,1+bν+m−(h),(c)+bν+2​m,(j)+bν+m;(−1)B+H+D+Kz]+∑ν=1HzhνΓ[(a)+hν+m,(g)+hν,(b)−hν+m,(h)ν′−hν(c)+hν+m,(j)+hν,(d)−hν+m,(k)−hν]×FB+H+C+J−1A+G+D+K[(a)+hν+m,(g)+hν,1−(d)+hν−m,1−(k)+hν;1−(b)+hν−m,1−(h)ν′+hν,(c)+hν+m,(j)+hν;(−1)B+H+D+Kz]}.\Sigma_{B}(z)=\sum_{m=0}^{\infty}\Gamma\Bigg{[}\begin{matrix}(e)+m,&(f)\\ (e),&(f)+m\end{matrix}\Bigg{]}\frac{x^{m}}{m!}\\ \times\Bigg{\{}\sum_{\nu=1}^{B}z^{b_{\nu}+m}\Gamma\Bigg{[}\begin{matrix}(a)+b_{\nu}+2m,&(g)+b_{\nu}+m,&(b)_{\nu}^{\prime}-b_{\nu},&(h)-b_{\nu}-m\\ (c)+b_{\nu}+2m,&(j)+b_{\nu}+m,&(d)-b_{\nu},&(k)-b_{\nu}-m\end{matrix}\Bigg{]}\times\\ {}_{A+G+D+K}F_{B+H+C+J-1}\Bigg{[}\begin{matrix}(a)+b_{\nu}+2m,&(g)+b_{\nu}+m,&1+b_{\nu}-(d),&1+b_{\nu}+m-(k);\\ 1+b_{\nu}-(b)_{\nu}^{\prime},&1+b_{\nu}+m-(h),&(c)+b_{\nu}+2m,&(j)+b_{\nu}+m;\end{matrix}\quad\\ (-1)^{B+H+D+K}z\Bigg{]}\\ +\sum_{\nu=1}^{H}z^{h_{\nu}}\Gamma\Bigg{[}\begin{matrix}(a)+h_{\nu}+m,&(g)+h_{\nu},&(b)-h_{\nu}+m,&(h)_{\nu}^{\prime}-h_{\nu}\\ (c)+h_{\nu}+m,&(j)+h_{\nu},&(d)-h_{\nu}+m,&(k)-h_{\nu}\end{matrix}\Bigg{]}\times\\ {}_{A+G+D+K}F_{B+H+C+J-1}\Bigg{[}\begin{matrix}(a)+h_{\nu}+m,&(g)+h_{\nu},&1-(d)+h_{\nu}-m,&1-(k)+h_{\nu};\\ 1-(b)+h_{\nu}-m,&1-(h)_{\nu}^{\prime}+h_{\nu},&(c)+h_{\nu}+m,&(j)+h_{\nu};\end{matrix}\\ \quad(-1)^{B+H+D+K}z\Bigg{]}\Bigg{\}}\,. (56)

These results can be transformed into Phchhammer symbols with a more elegant form,

ΣA(z)=∑μ=1AΓ[(a)μ′−aμ,(b)+aμ,(g)−aμ,(h)+aμ(c)−aμ,(d)+aμ,(j)−aμ,(k)+aμ]×∑m=0∞∑n=0∞((b)+aμ)2​m+n​((h)+aμ)m+n​((e))m​(1+aμ−(c))n​(1+aμ−(j))m+n(1+aμ−(a)μ′)n​(1+aμ−(g))m+n​((f))m​((d)+aμ)2​m+n​((k)+aμ)m+n⋅xm​z−aμ−m−n​(−1)n​(A+G−C−J)m!​n!×Γ​((j)−aμ)​Γ​(1+aμ−(j))Γ​(1+aμ−(j)+m)​Γ​((j)−aμ−m)⋅Γ​(1+aμ−(g)+m)​Γ​((g)−aμ−m)Γ​((g)−aμ)​Γ​(1+aμ−(g))+∑μ=1GΓ[(a)−gμ,(b)+gμ,(g)μ′−gμ,(h)+gμ(c)−gμ,(d)+gμ,(j)−gμ,(k)+gμ]×∑m=0∞∑n=0∞((a)−gμ)m−n​((b)+gμ)m+n​(1+gμ−(j))n​((h)+gμ)n​((e))m((c)−gμ)m−n​((d)+gμ)m+n​(1+gμ−(g)μ′)n​((k)+gμ)n​((f))m⋅xm​z−gμ−n​(−1)n​(A+G−C−J)m!​n!×Γ​((a)−gμ+m)​Γ​(1+gμ−(a)−m)Γ​((a)−gμ+m−n)​Γ​(1+gμ−(a)−m+n)×Γ​((c)−gμ+m−n)​Γ​(1+gμ−(c)−m+n)Γ​((c)−gμ+m)​Γ​(1+gμ−(c)−m),\Sigma_{A}(z)=\sum_{\mu=1}^{A}\Gamma\Bigg{[}\begin{matrix}(a)_{\mu}^{\prime}-a_{\mu},&(b)+a_{\mu},&(g)-a_{\mu},&(h)+a_{\mu}\\ (c)-a_{\mu},&(d)+a_{\mu},&(j)-a_{\mu},&(k)+a_{\mu}\end{matrix}\Bigg{]}\times\\ \sum_{m=0}^{\infty}\sum_{n=0}^{\infty}\frac{((b)+a_{\mu})_{2m+n}((h)+a_{\mu})_{m+n}((e))_{m}(1+a_{\mu}-(c))_{n}(1+a_{\mu}-(j))_{m+n}}{(1+a_{\mu}-(a)_{\mu}^{\prime})_{n}(1+a_{\mu}-(g))_{m+n}((f))_{m}((d)+a_{\mu})_{2m+n}((k)+a_{\mu})_{m+n}}\\ \cdot\frac{x^{m}z^{-a_{\mu}-m-n}(-1)^{n(A+G-C-J)}}{m!n!}\\ \times\frac{\Gamma((j)-a_{\mu})\Gamma(1+a_{\mu}-(j))}{\Gamma(1+a_{\mu}-(j)+m)\Gamma((j)-a_{\mu}-m)}\cdot\frac{\Gamma(1+a_{\mu}-(g)+m)\Gamma((g)-a_{\mu}-m)}{\Gamma((g)-a_{\mu})\Gamma(1+a_{\mu}-(g))}\\ +\sum_{\mu=1}^{G}\Gamma\Bigg{[}\begin{matrix}(a)-g_{\mu},&(b)+g_{\mu},&(g)_{\mu}^{\prime}-g_{\mu},&(h)+g_{\mu}\\ (c)-g_{\mu},&(d)+g_{\mu},&(j)-g_{\mu},&(k)+g_{\mu}\end{matrix}\Bigg{]}\times\\ \sum_{m=0}^{\infty}\sum_{n=0}^{\infty}\frac{((a)-g_{\mu})_{m-n}((b)+g_{\mu})_{m+n}(1+g_{\mu}-(j))_{n}((h)+g_{\mu})_{n}((e))_{m}}{((c)-g_{\mu})_{m-n}((d)+g_{\mu})_{m+n}(1+g_{\mu}-(g)_{\mu}^{\prime})_{n}((k)+g_{\mu})_{n}((f))_{m}}\\ \cdot\frac{x^{m}z^{-g_{\mu}-n}(-1)^{n(A+G-C-J)}}{m!n!}\\ \times\frac{\Gamma((a)-g_{\mu}+m)\Gamma(1+g_{\mu}-(a)-m)}{\Gamma((a)-g_{\mu}+m-n)\Gamma(1+g_{\mu}-(a)-m+n)}\\ \times\frac{\Gamma((c)-g_{\mu}+m-n)\Gamma(1+g_{\mu}-(c)-m+n)}{\Gamma((c)-g_{\mu}+m)\Gamma(1+g_{\mu}-(c)-m)}, (57)

and

ΣB(z)=∑ν=1BΓ[(a)+bν,(b)ν′−bν,(g)+bν,(h)−bν(c)+bν,(d)−bν,(j)+bν,(k)−bν]×∑m=0∞∑n=0∞((a)+bν)2​m+n​((g)+bν)m+n​((e))m​(1+bν−(d))n​(1+bν−(k))m+n(1+bν−(b)ν′)n​(1+bν−(h))m+n​((f))m​((c)+bν)2​m+n​((j)+bν)m+n⋅xm​zbν+m+n​(−1)n​(B+H−K−D)m!​n!×Γ​((h)−bν−m)​Γ​(1+bν+m−(h))Γ​(1+bν−(h))​Γ​((h)−bν)⋅Γ​(1+bν−(k))​Γ​((k)−bν)Γ​((k)−bν−m)​Γ​(1+bν+m−(k))+∑ν=1HΓ[(a)+hν,(b)−hν,(g)+hν,(h)ν′−hν(c)+hν,(d)−hν,(j)+hν,(k)−hν]×∑m=0∞∑n=0∞((a)+hν)m+n​((b)−hν)m−n​(1+hν−(k))n​((g)+hν)n​((e))m((d)−hν)m−n​((c)+hν)m+n​(1+hν−(h)ν′)n​((j)+hν)n​((f))m⋅xm​zhν+n​(−1)n​(B+H−K−D)m!​n!×Γ​((b)−hν+m)​Γ​(1−(b)+hν−m)Γ​((b)−hν+m−n)​Γ​(1−(b)+hν−m+n)×Γ​((d)−hν+m−n)​Γ​(1−(d)+hν−m+n)Γ​((d)−hν+m)​Γ​(1−(d)+hν−m).\Sigma_{B}(z)=\sum_{\nu=1}^{B}\Gamma\Bigg{[}\begin{matrix}(a)+b_{\nu},&(b)_{\nu}^{\prime}-b_{\nu},&(g)+b_{\nu},&(h)-b_{\nu}\\ (c)+b_{\nu},&(d)-b_{\nu},&(j)+b_{\nu},&(k)-b_{\nu}\end{matrix}\Bigg{]}\times\\ \sum_{m=0}^{\infty}\sum_{n=0}^{\infty}\frac{((a)+b_{\nu})_{2m+n}((g)+b_{\nu})_{m+n}((e))_{m}(1+b_{\nu}-(d))_{n}(1+b_{\nu}-(k))_{m+n}}{(1+b_{\nu}-(b)_{\nu}^{\prime})_{n}(1+b_{\nu}-(h))_{m+n}((f))_{m}((c)+b_{\nu})_{2m+n}((j)+b_{\nu})_{m+n}}\\ \cdot\frac{x^{m}z^{b_{\nu}+m+n}(-1)^{n(B+H-K-D)}}{m!n!}\\ \times\frac{\Gamma((h)-b_{\nu}-m)\Gamma(1+b_{\nu}+m-(h))}{\Gamma(1+b_{\nu}-(h))\Gamma((h)-b_{\nu})}\cdot\frac{\Gamma(1+b_{\nu}-(k))\Gamma((k)-b_{\nu})}{\Gamma((k)-b_{\nu}-m)\Gamma(1+b_{\nu}+m-(k))}\\ +\sum_{\nu=1}^{H}\Gamma\Bigg{[}\begin{matrix}(a)+h_{\nu},&(b)-h_{\nu},&(g)+h_{\nu},&(h)_{\nu}^{\prime}-h_{\nu}\\ (c)+h_{\nu},&(d)-h_{\nu},&(j)+h_{\nu},&(k)-h_{\nu}\end{matrix}\Bigg{]}\times\\ \sum_{m=0}^{\infty}\sum_{n=0}^{\infty}\frac{((a)+h_{\nu})_{m+n}((b)-h_{\nu})_{m-n}(1+h_{\nu}-(k))_{n}((g)+h_{\nu})_{n}((e))_{m}}{((d)-h_{\nu})_{m-n}((c)+h_{\nu})_{m+n}(1+h_{\nu}-(h)_{\nu}^{\prime})_{n}((j)+h_{\nu})_{n}((f))_{m}}\\ \cdot\frac{x^{m}z^{h_{\nu}+n}(-1)^{n(B+H-K-D)}}{m!n!}\\ \times\frac{\Gamma((b)-h_{\nu}+m)\Gamma(1-(b)+h_{\nu}-m)}{\Gamma((b)-h_{\nu}+m-n)\Gamma(1-(b)+h_{\nu}-m+n)}\\ \times\frac{\Gamma((d)-h_{\nu}+m-n)\Gamma(1-(d)+h_{\nu}-m+n)}{\Gamma((d)-h_{\nu}+m)\Gamma(1-(d)+h_{\nu}-m)}. (58)

To simplify the redundant Γ\Gamma functions appearing in the last line of each summation series, we make use of Euler’s reflection formula on Γ\Gamma functions.

Lemma A.1 (Euler’s reflection formula).
Γ​(1−z)​Γ​(z)=πsin⁡π​z,z∉ℤ.\Gamma(1-z)\Gamma(z)=\frac{\pi}{\sin\pi z},\qquad z\notin\mathbb{Z}. (59)

A concise proof delivered by Dedekind can be found in ref. srinivasan2011dedekind . Replacing zz by z−nz-n in Euler’s reflection formula and noting that

πsin⁡π​(z−n)=(−1)n​πsin⁡π​z,\frac{\pi}{\sin\pi(z-n)}=(-1)^{n}\frac{\pi}{\sin\pi z}, (60)

we derived another useful lemma.

Lemma A.2.

For n∈ℤn\in\mathbb{Z} and z∉ℤz\notin\mathbb{Z},

Γ​(1−z)​Γ​(z)Γ​(1+n−z)​Γ​(z−n)=(−1)n.\frac{\Gamma(1-z)\Gamma(z)}{\Gamma(1+n-z)\Gamma(z-n)}=(-1)^{n}\,. (61)

In our case, the redundant Γ\Gamma functions cancel out,

Γ​((j)−aμ)​Γ​(1+aμ−(j))Γ​(1+aμ−(j)+m)​Γ​((j)−aμ−m)\displaystyle\frac{\Gamma((j)-a_{\mu})\Gamma(1+a_{\mu}-(j))}{\Gamma(1+a_{\mu}-(j)+m)\Gamma((j)-a_{\mu}-m)} =(−1)m​J,\displaystyle=(-1)^{mJ}, (62)
Γ​(1+aμ−(g)+m)​Γ​((g)−aμ−m)Γ​((g)−aμ)​Γ​(1+aμ−(g))\displaystyle\frac{\Gamma(1+a_{\mu}-(g)+m)\Gamma((g)-a_{\mu}-m)}{\Gamma((g)-a_{\mu})\Gamma(1+a_{\mu}-(g))} =(−1)m​G,\displaystyle=(-1)^{mG}, (63)
Γ​((a)−gμ+m)​Γ​(1+gμ−(a)−m)Γ​((a)−gμ+m−n)​Γ​(1+gμ−(a)−m+n)\displaystyle\frac{\Gamma((a)-g_{\mu}+m)\Gamma(1+g_{\mu}-(a)-m)}{\Gamma((a)-g_{\mu}+m-n)\Gamma(1+g_{\mu}-(a)-m+n)} =(−1)n​A,\displaystyle=(-1)^{nA}, (64)
Γ​((c)−gμ+m−n)​Γ​(1+gμ−(c)−m+n)Γ​((c)−gμ+m)​Γ​(1+gμ−(c)−m)\displaystyle\frac{\Gamma((c)-g_{\mu}+m-n)\Gamma(1+g_{\mu}-(c)-m+n)}{\Gamma((c)-g_{\mu}+m)\Gamma(1+g_{\mu}-(c)-m)} =(−1)n​C,\displaystyle=(-1)^{nC}, (65)
Γ​((h)−bν−m)​Γ​(1+bν+m−(h))Γ​(1+bν−(h))​Γ​((h)−bν)\displaystyle\frac{\Gamma((h)-b_{\nu}-m)\Gamma(1+b_{\nu}+m-(h))}{\Gamma(1+b_{\nu}-(h))\Gamma((h)-b_{\nu})} =(−1)m​H,\displaystyle=(-1)^{mH}, (66)
Γ​(1+bν−(k))​Γ​((k)−bν)Γ​((k)−bν−m)​Γ​(1+bν+m−(k))\displaystyle\frac{\Gamma(1+b_{\nu}-(k))\Gamma((k)-b_{\nu})}{\Gamma((k)-b_{\nu}-m)\Gamma(1+b_{\nu}+m-(k))} =(−1)m​K,\displaystyle=(-1)^{mK}, (67)
Γ​((b)−hν+m)​Γ​(1−(b)+hν−m)Γ​((b)−hν+m−n)​Γ​(1−(b)+hν−m+n)\displaystyle\frac{\Gamma((b)-h_{\nu}+m)\Gamma(1-(b)+h_{\nu}-m)}{\Gamma((b)-h_{\nu}+m-n)\Gamma(1-(b)+h_{\nu}-m+n)} =(−1)n​B,\displaystyle=(-1)^{nB}, (68)
Γ​((d)−hν+m−n)​Γ​(1−(d)+hν−m+n)Γ​((d)−hν+m)​Γ​(1−(d)+hν−m)\displaystyle\frac{\Gamma((d)-h_{\nu}+m-n)\Gamma(1-(d)+h_{\nu}-m+n)}{\Gamma((d)-h_{\nu}+m)\Gamma(1-(d)+h_{\nu}-m)} =(−1)n​D.\displaystyle=(-1)^{nD}. (69)

Therefore, we obtain

ΣA(z)=∑μ=1AΓ[(a)μ′−aμ,(b)+aμ,(g)−aμ,(h)+aμ(c)−aμ,(d)+aμ,(j)−aμ,(k)+aμ]×∑m=0∞∑n=0∞((b)+aμ)2​m+n​((h)+aμ)m+n​((e))m​(1+aμ−(c))n​(1+aμ−(j))m+n(1+aμ−(a)μ′)n​(1+aμ−(g))m+n​((f))m​((d)+aμ)2​m+n​((k)+aμ)m+n⋅xm​z−aμ−m−n​(−1)n​(A+G−C−J)+m​(G−J)m!​n!+∑μ=1GΓ[(a)−gμ,(b)+gμ,(g)μ′−gμ,(h)+gμ(c)−gμ,(d)+gμ,(j)−gμ,(k)+gμ]×∑m=0∞∑n=0∞((a)−gμ)m−n​((b)+gμ)m+n​(1+gμ−(j))n​((h)+gμ)n​((e))m((c)−gμ)m−n​((d)+gμ)m+n​(1+gμ−(g)μ′)n​((k)+gμ)n​((f))m⋅xm​z−gμ−n​(−1)n​(G−J)m!​n!,\Sigma_{A}(z)=\sum_{\mu=1}^{A}\Gamma\Bigg{[}\begin{matrix}(a)_{\mu}^{\prime}-a_{\mu},&(b)+a_{\mu},&(g)-a_{\mu},&(h)+a_{\mu}\\ (c)-a_{\mu},&(d)+a_{\mu},&(j)-a_{\mu},&(k)+a_{\mu}\end{matrix}\Bigg{]}\times\\ \sum_{m=0}^{\infty}\sum_{n=0}^{\infty}\frac{((b)+a_{\mu})_{2m+n}((h)+a_{\mu})_{m+n}((e))_{m}(1+a_{\mu}-(c))_{n}(1+a_{\mu}-(j))_{m+n}}{(1+a_{\mu}-(a)_{\mu}^{\prime})_{n}(1+a_{\mu}-(g))_{m+n}((f))_{m}((d)+a_{\mu})_{2m+n}((k)+a_{\mu})_{m+n}}\\ \cdot\frac{x^{m}z^{-a_{\mu}-m-n}(-1)^{n(A+G-C-J)+m(G-J)}}{m!n!}\\ +\sum_{\mu=1}^{G}\Gamma\Bigg{[}\begin{matrix}(a)-g_{\mu},&(b)+g_{\mu},&(g)_{\mu}^{\prime}-g_{\mu},&(h)+g_{\mu}\\ (c)-g_{\mu},&(d)+g_{\mu},&(j)-g_{\mu},&(k)+g_{\mu}\end{matrix}\Bigg{]}\times\\ \sum_{m=0}^{\infty}\sum_{n=0}^{\infty}\frac{((a)-g_{\mu})_{m-n}((b)+g_{\mu})_{m+n}(1+g_{\mu}-(j))_{n}((h)+g_{\mu})_{n}((e))_{m}}{((c)-g_{\mu})_{m-n}((d)+g_{\mu})_{m+n}(1+g_{\mu}-(g)_{\mu}^{\prime})_{n}((k)+g_{\mu})_{n}((f))_{m}}\\ \cdot\frac{x^{m}z^{-g_{\mu}-n}(-1)^{n(G-J)}}{m!n!}, (70)

and

ΣB(z)=∑ν=1BΓ[(a)+bν,(b)ν′−bν,(g)+bν,(h)−bν(c)+bν,(d)−bν,(j)+bν,(k)−bν]×∑m=0∞∑n=0∞((a)+bν)2​m+n​((g)+bν)m+n​((e))m​(1+bν−(d))n​(1+bν−(k))m+n(1+bν−(b)ν′)n​(1+bν−(h))m+n​((f))m​((c)+bν)2​m+n​((j)+bν)m+n⋅xm​zbν+m+n​(−1)n​(B+H−K−D)+m​(H−K)m!​n!+∑ν=1HΓ[(a)+hν,(b)−hν,(g)+hν,(h)ν′−hν(c)+hν,(d)−hν,(j)+hν,(k)−hν]×∑m=0∞∑n=0∞((a)+hν)m+n​((b)−hν)m−n​(1+hν−(k))n​((g)+hν)n​((e))m((d)−hν)m−n​((c)+hν)m+n​(1+hν−(h)ν′)n​((j)+hν)n​((f))m⋅xm​zhν+n​(−1)n​(H−K)m!​n!.\Sigma_{B}(z)=\sum_{\nu=1}^{B}\Gamma\Bigg{[}\begin{matrix}(a)+b_{\nu},&(b)_{\nu}^{\prime}-b_{\nu},&(g)+b_{\nu},&(h)-b_{\nu}\\ (c)+b_{\nu},&(d)-b_{\nu},&(j)+b_{\nu},&(k)-b_{\nu}\end{matrix}\Bigg{]}\times\\ \sum_{m=0}^{\infty}\sum_{n=0}^{\infty}\frac{((a)+b_{\nu})_{2m+n}((g)+b_{\nu})_{m+n}((e))_{m}(1+b_{\nu}-(d))_{n}(1+b_{\nu}-(k))_{m+n}}{(1+b_{\nu}-(b)_{\nu}^{\prime})_{n}(1+b_{\nu}-(h))_{m+n}((f))_{m}((c)+b_{\nu})_{2m+n}((j)+b_{\nu})_{m+n}}\\ \cdot\frac{x^{m}z^{b_{\nu}+m+n}(-1)^{n(B+H-K-D)+m(H-K)}}{m!n!}\\ +\sum_{\nu=1}^{H}\Gamma\Bigg{[}\begin{matrix}(a)+h_{\nu},&(b)-h_{\nu},&(g)+h_{\nu},&(h)_{\nu}^{\prime}-h_{\nu}\\ (c)+h_{\nu},&(d)-h_{\nu},&(j)+h_{\nu},&(k)-h_{\nu}\end{matrix}\Bigg{]}\times\\ \sum_{m=0}^{\infty}\sum_{n=0}^{\infty}\frac{((a)+h_{\nu})_{m+n}((b)-h_{\nu})_{m-n}(1+h_{\nu}-(k))_{n}((g)+h_{\nu})_{n}((e))_{m}}{((d)-h_{\nu})_{m-n}((c)+h_{\nu})_{m+n}(1+h_{\nu}-(h)_{\nu}^{\prime})_{n}((j)+h_{\nu})_{n}((f))_{m}}\\ \cdot\frac{x^{m}z^{h_{\nu}+n}(-1)^{n(H-K)}}{m!n!}. (71)

which are the same results as in Corollary 2.1.

Appendix B The explicit expressions of 𝒜i,ℬi\mathcal{A}_{i},\mathcal{B}_{i} and 𝒞i\mathcal{C}_{i} in eq. (43).

𝒜1\displaystyle\mathcal{A}_{1} =∑i=215σi+D2,\displaystyle=\sum_{i=2}^{15}\sigma_{i}+\frac{D}{2}, (72)
𝒜2\displaystyle\mathcal{A}_{2} =σ2,4,5,6,9,11,14,16+ν1,\displaystyle=\sigma_{2,4,5,6,9,11,14,16}+\nu_{1}, (73)
𝒜3\displaystyle\mathcal{A}_{3} =σ3,4,7,9,10,12,15,17+ν2,\displaystyle=\sigma_{3,4,7,9,10,12,15,17}+\nu_{2}, (74)
𝒜4\displaystyle\mathcal{A}_{4} =1+σ2,3,9,10,11,12,13−σ20+D2−ν5,\displaystyle=1+\sigma_{2,3,9,10,11,12,13}-\sigma_{20}+\frac{D}{2}-\nu_{5}, (75)
𝒜5\displaystyle\mathcal{A}_{5} =1+σ2,3,4,5,6,7,8,16,17,18,19,20+ν−D−ν6,\displaystyle=1+\sigma_{2,3,4,5,6,7,8,16,17,18,19,20}+\nu-D-\nu_{6}, (76)
ℬ1\displaystyle\mathcal{B}_{1} =σ4,5,6,7,8,14,15,20+ν5,\displaystyle=\sigma_{4,5,6,7,8,14,15,20}+\nu_{5}, (77)
ℬ2\displaystyle\mathcal{B}_{2} =−σ3,10,12,13+σ4,5,6,14,16,20−D2+ν1+ν5,\displaystyle=-\sigma_{3,10,12,13}+\sigma_{4,5,6,14,16,20}-\frac{D}{2}+\nu_{1}+\nu_{5}, (78)
ℬ3\displaystyle\mathcal{B}_{3} =−σ2,11,13+σ4,7,15,17,20−D2+ν2+ν5,\displaystyle=-\sigma_{2,11,13}+\sigma_{4,7,15,17,20}-\frac{D}{2}+\nu_{2}+\nu_{5}, (79)
ℬ4\displaystyle\mathcal{B}_{4} =1−σ2,3,9,10,11,12,13+σ20−D2+ν5,\displaystyle=1-\sigma_{2,3,9,10,11,12,13}+\sigma_{20}-\frac{D}{2}+\nu_{5}, (80)
ℬ5\displaystyle\mathcal{B}_{5} =1−σ4,5,6,7,8,16,17,19+σ9,10,11,12,13−2​σ20−ν+3​D2+ν6−ν5,\displaystyle=1-\sigma_{4,5,6,7,8,16,17,19}+\sigma_{9,10,11,12,13}-2\sigma_{20}-\nu+\frac{3D}{2}+\nu_{6}-\nu_{5}, (81)
𝒞1\displaystyle\mathcal{C}_{1} =σ9,10,11,12,13,14,15−σ16,17,18,19,20−ν+3​D2+ν6,\displaystyle=\sigma_{9,10,11,12,13,14,15}-\sigma_{16,17,18,19,20}-\nu+\frac{3D}{2}+\nu_{6}, (82)
𝒞2\displaystyle\mathcal{C}_{2} =−σ3,7,8,17,18,19,20+σ9,11,14−ν+D+ν1+ν6,\displaystyle=-\sigma_{3,7,8,17,18,19,20}+\sigma_{9,11,14}-\nu+D+\nu_{1}+\nu_{6}, (83)
𝒞3\displaystyle\mathcal{C}_{3} =−σ2,5,6,8,16,18,19,20−+σ9,10,12,15−ν+D+ν2+ν6,\displaystyle=-\sigma_{2,5,6,8,16,18,19,20}-+\sigma_{9,10,12,15}-\nu+D+\nu_{2}+\nu_{6}, (84)
𝒞4\displaystyle\mathcal{C}_{4} =1+σ2−σ3,4,5,6,7,8,16,17,18,19,20−ν+D+ν6,\displaystyle=1+\sigma_{2}-\sigma_{3,4,5,6,7,8,16,17,18,19,20}-\nu+D+\nu_{6}, (85)
𝒞5\displaystyle\mathcal{C}_{5} =1−σ4,5,6,7,8,16,17,18,19+σ9,10,11,12,13−2​σ20−ν−3​D2−ν5+ν6.\displaystyle=1-\sigma_{4,5,6,7,8,16,17,18,19}+\sigma_{9,10,11,12,13}-2\sigma_{20}-\nu-\frac{3D}{2}-\nu_{5}+\nu_{6}. (86)

References