Hypergeometric Series Representations of Feynman Integrals by GKZ Hypergeometric Systems

We show that almost all Feynman integrals as well as their coefficients in a Laurent series in dimensional regularization can be written in terms of Horn hypergeometric functions. By applying the results of Gelfand-Kapranov-Zelevinsky (GKZ) we derive a formula for a class of hypergeometric series representations of Feynman integrals, which can be obtained by triangulations of the Newton polytope $\Delta_G$ corresponding to the Lee-Pomeransky polynomial $G$. Those series can be of higher dimension, but converge fast for convenient kinematics, which also allows numerical applications. Further, we discuss possible difficulties which can arise in a practical usage of this approach and give strategies to solve them.


Introduction
In the early days of calculating Feynman amplitudes, it was proposed by Regge to consider Feynman integrals as a kind of generalized hypergeometric functions [60], where the singularities of those hypergeometric functions coincide with the Landau singularities. Later on Kashiwara and Kawai [43] showed that Feynman integrals satisfy indeed holonomic differential equations, where the singularities of those holonomic differential equations are determined by the Landau singularities.
Apart from characterizing the Feynman integral by "hypergeometric" partial differential equation systems, many applications determine the Feynman integral as a generalized hypergeometric series. Usually, the often used Mellin-Barnes approach [68] results in Pochhammer series p F q , Appell functions, Lauricella functions and related functions by applying the residue theorem [8]. Furthermore, for arbitrary one-loop Feynman integrals it is known that they can always be represented by a small set of hypergeometric series [24]. Thirdly, the Feynman integral may be expressed by "hypergeometric" integrals like the generalized Meijer G-or Fox H-functions [15,39,40].
Thus there arise three different notions of the term "hypergeometric" in the Feynman integral calculus, where every notion generalizes different characterizations of the classical hypergeometric Gauß function 2 F 1 (a, b, c; x). In the late 1980s Gelfand, Kapranov, Zelevinsky (GKZ) and collaborators [26][27][28][29][30][31] were starting to develop a comprehensive method to generalize the notion of "hypergeometric" functions in a consistent way. Those functions are called A-hypergeometric functions and are defined by a special holonomic system of partial differential equations.
As Gelfand, Kapranov and Zelevinsky illustrated with Euler integrals, the GKZ approach not only generalizes the concept of hypergeometric functions but can also be used for analyzing and solving integrals [30].
For physicists the GKZ perspective is not entirely new. Already in the 1990s, string theorists applied the GKZ approach in order to calculate period integrals and worked out the mirror symmetry [37,38]. Recently, the GKZ approach was also used to obtain differential equations for the Feynman integral from the maximal cut [74]. Still, the approach of Gelfand, Kapranov and Zelevinsky is no common practice among physicists.
In 2016 Nasrollahpoursamami showed that the Feynman integral satisfies a differential equation system which is isomorphic to a GKZ system [52]. Very recently 1 , this fact was also shown directly in [20] based on the Lee-Pomeransky representation of the Feynman integral.
Beyond the above statement which characterizes generic Feynman integrals as Ahypergeometric functions, we show that generic Feynman integrals as well as every coefficient in the expansion in dimensional regularization belong to the class of Horn hypergeometric functions. Furthermore, we give an explicit formula for a multivariate series representation of a generic Feynman integral for unimodular triangulations. This allows, to evaluate the Feynman integral efficiently for convenient kinematic regions.
Therefore, from the perspective of GKZ it turns out that Horn hypergeometric functions and Feynman integrals share many properties, e.g. Horn hypergeometric functions satisfy special relations similar to the IBP-relations of the Feynman integral [47]. In this article, we work out the connection between Feynman integrals and Horn hypergeometric functions from the GKZ perspective and we give a strategy, as to how to use this knowledge in the evaluation of Feynman integrals.
The connections between Feynman integrals and hypergeometric functions was investigated over decades and a comprehensive summary of these investigations can be found in [41]. Horn hypergeometric functions also appear often in the Mellin-Barnes approach and have been studied intensively by other authors, e.g [16,61].
The paper is structured as follows. After recalling the parametric representations of the Feynman integral, which are crucial for this approach, we mention some properties of the Feynman integral, which can be obtained from their algebraic geometry description. This first section will also include a short mathematical interlude about convex polytopes, which are necessary in the description of the properties of Feynman integrals, as well as in the later GKZ approach. Secondly, we briefly introduce the GKZ hypergeometric system, which is a system of partial differential equations and recall a series solution of those systems. This leads us to the last section, in which we merge these two aspects in order to derive an analytical series representation of the Feynman integral. To complete this section we also discuss some features as well as possible difficulties which can arise in the evaluation of Feynman integrals with GKZ systems. We conclude this article by calculating the full massive sunset Feynman diagram as a non-trivial example in order to illustrate this procedure and to give a glimpse of its further scope of application.
Quantum Field Theory. I would like to express my special thanks to Christian Bogner for his great encouragement, his supervision and his assistance in difficult circumstances. Further, I would like to thank Dirk Kreimer for his support, organizational skills and for including me in his group. I am also grateful to Josua Faller and Tatiana Stahlhut for their proofreading and all the members of Dirk Kreimer's group for helpful discussions.

Feynman Integrals
In this section we shortly recall the parametric representations of Feynman integrals, which are based on graph polynomials [7]. Those parametric representations contain normally two graph polynomials, known as the Symanzik polynomials. Recently, Lee and Pomeransky [48] introduced a parametric representation, which depends only on one graph polynomial, which simplifies the application of the Gelfand-Kapranov-Zelevinsky approach. In this representation the Feynman integral can be formally described as a multivariate Euler-Mellin integral.
In the second part of this section, we rely on some general properties of Euler-Mellin integrals in order to recall some basic properties of the Feynman integral from the perspective of the Lee-Pomeransky representation. To keep this discussion short, we refer to standard literature [68,70,77] for further properties of the Feynman integral.
Since the formulation of this properties of Feynman integrals, as well as the later introduced GKZ approach, include convex polytopes, we give a short mathematical interlude about convex polytopes between these two parts. A proceeding and more detailed description of polytopes can be found in the appendix.

Parametric Representations of Feynman Integrals
The Feynman integral in momentum space is an (Ld)-dimensional loop integral over propagators. Since propagators are at most quadratic in the loop momenta, one can rewrite the Feynman integral as an n-dimensional integral over Schwinger parameters x 1 , . . . , x n ∈ R ≥0 [57]. With this rephrasing one can make the Feynman integral also meaningful for d ∈ C, as required in the procedure of dimensional regularization [73].
In the parametric version of Feynman integrals, two polynomials in the Schwinger parameters arise, which are known as first and second Symanzik polynomial [7] U = Here, T i is the set of spanning forests of the Feynman graph Γ consisting of i components, E denotes the set of all edges of Γ and s F is the squared momentum flowing from the one component of the 2-forest F to the other component. For Euclidean kinematics the coefficients of all monomials in the second Symanzik polynomial will be positive. In the following we restrict the discussion to Euclidean kinematics in order to avoid zeros of the Symanzik polynomials in the positive orthant x ∈ R n >0 . Since almost every Feynman integral with Minkowskian kinematics has a non-vanishing overlap with the Euclidean region, one can extend the Euclidean result to the Minkowskian result by analytic continuation. In this procedure additional divergences can appear and thus the analytic continuation can be far from trivial [57].
As the starting point of this discussion we define the Feynman integral in the Feynman parametric representation. Definition 1.1 [Feynman integral]: For a given Feynman graph Γ with n edges, L loops 2 and the Symanzik polynomials U and F , the Feynman integral is given as an n-dimensional integral where ω := n i=1 ν i − Ld 2 is the superficial degree of divergence and δ(x) denotes the Dirac δ function. The ν ∈ C n with Re ν i > 0 are the propagator powers. Note that for simplicity we use a multi-index notation throughout the whole paper, which implies the following shorthand notations Clearly, those integrals converge not for every values of ν ∈ C n and d ∈ C and also not for every choice of the polynomials U and F . It can be shown, that besides of the class of massless tadpole graphs, these integrals define meromorphic functions in ν and d. Therefore, we will consider the Feynman integral always as the meromorphic continuation of the integral (1.3) to the whole complex plane. The convergence as well as the meromorphic expansion will be discussed in more detail in section 1.3.

Remark:
One is often interested in ν i ∈ N. However, in the following it will be convenient to consider a slightly more general notion of Feynman integrals where the propagator powers ν i are not restricted to integer values only. Just the restriction Re ν i > 0 is necessary to guaranty the parametric rewriting and convergence. These conditions are often assumed in the literature e.g. in [6,69].
The parametric representation in equation (1.3) is not the only representation in terms of Schwinger parameters. For the following approach another parametric representation invented by Lee and Pomeransky [48] is more convenient. [48]]: The Feynman integral from equation (1.3) can also be written as where the integral depends only on the sum of Symanzik polynomials G = U + F , which we call the Lee-Pomeransky polynomial G. The equality of representations is in the sense of meromorphic extension.
Proof. (A proof can be found also in [6]) Since U is homogeneous of degree L and F is homogeneous of degree L + 1, the integral in (1.5) as a function of d converges in the strip As we consider Re ν i > 0 it is Λ = ∅. For an equality in the sense of meromorphic extension it is sufficient to show that there is a non-vanishing interval where the equality of representations holds. (1.5), changing the integration order and substituting x i → sx i one obtains (1.6) can be calculated explicitly as a beta function in the region Λ (and for U, F > 0) one attains the representation (1.3).
Thus, the whole structure of a Feynman integral can be expressed in only one single polynomial G ∈ C[x 1 , . . . , x n ] . As it is possible for every polynomial, one can write the Lee-Pomeransky polynomial G as where A is a finite set consisting of N pairwise distinct column vectors a j ∈ Z n ≥0 and the coefficients of this polynomial are complex numbers z ∈ (C \ {0}) N , which contain the kinematics and masses. This representation is unique, up to the possibility of different monomial orderings. Without loss of generality we can fix an arbitrary monomial ordering and we denote the set of column vectors A in a matrix structure. As we will observe later, the exponential vectors of the Lee-Pomeransky polynomial G z satisfy an affine structure. Thus, it will be convenient to define Further, for an index set σ ⊂ {1, . . . , N } the notation A σ denotes the restriction of A according to columns indexed by σ.
In addition, we define coordinates ν := (ν 0 , ν) ∈ C n+1 , where the first entry of ν contains the physical dimension ν 0 := d 2 . From a mathematical point of view, the physical dimension d and the propagator powers ν i have a similar role in parametric Feynman integrals. This is the reason why instead of using dimensional regularization one can also regularize the integral by considering the propagator powers as complex numbers as done in analytic regularization [69].
In fact, this definition of G z (x) includes a generalization of the original Feynman integral, since also the first Symanzik polynomial gets coefficients z i . Thus, we define: Definition 1.2 [Generic Feynman Integrals]: Let G z (x) a Lee-Pomeransky polynomial with generic coefficients z ∈ C N satisfying Re z j > 0. The generic Feynman integral is the meromorphic continuation of the integral Remark: In this definition the complete graph structure, which is necessary to evaluate the Feynman integral, is given by the matrix A. The variables z ∈ C N contain the physical information about kinematics and masses and the ν ∈ C n+1 are the regularization parameters.
In contrast to equation (1.5) the coefficients in G z which come from the first Symanzik polynomial are treated as generic, instead of equal to 1. We will discuss later how one can remove these auxiliary variables afterwards to obtain the "ordinary" Feynman integral. To avoid unnecessary prefactors in the following, we omit also the factor Γ(ν 0 − ω)Γ(ν) in this definition With this definitions one can derive another representation of the Feynman integral as a multi-dimensional Mellin-Barnes integral.
Theorem 1.2 [Representation as Fox H-function]: Let σ ⊂ {1, . . . , N } be an index subset with cardinality n+1, such that the matrix A restricted to columns of σ is invertible, det A σ = 0. Then the Feynman integral can be written as the multi-dimensional Mellin-Barnes integral wherever this integral converges. The setσ := {1, . . . N } \ σ denotes the complement of σ, containing r := N − n − 1 elements. Restrictions of vectors and matrices to those index sets are similarly defined as z σ := (z i ) i∈σ , zσ := (z i ) i∈σ , Aσ := (a i ) i∈σ . Every component of the integration contour γ ∈ C r goes from −i∞ to i∞ such that the poles of the integrand are separated. Corollary 1.3: Let N = n + 1 or in other words let A be quadratic. If there is a region D ⊆ C n+1 such that the Feynman integral J A (ν, z) converges absolutely for ν ∈ D, the matrix A is invertible and the Feynman integral is only a simple combination of Γfunctions: Proof. Starting from equation (1.9) by the Schwinger trick one gets where ν 0 = d 2 . Writing x = (x 0 , x), ν = (ν 0 , ν) and using the Cahen-Mellin integral representation of exponential function one obtains with u ∈ C n+1 , t ∈ C r , some arbitrary positive numbers δ i > 0 and where we split the polynomial G into a σ and aσ part. By a substitution u → A −1 σ u it is (1.14) Since the matrix A σ contains only positive values, the integration region remains the same A σ δ + iA σ R n+1 δ + iR n+1 with some other positive numbers δ ∈ R n+1 >0 , which additionally have to satisfy A −1 σ δ > 0. By Mellin's inversion theorem [1] only the tintegration remains and one obtains equation (1.10).
Thereby, the integration contour has to be chosen, such that the poles are separated from each other in order to satisfy A −1 σ δ > 0. More specific this means that the contour γ has the form c + iR n+1 where c ∈ R n+1 >0 satisfies A −1 σ ν − A −1 σ Aσc > 0. Clearly, in order for those c to exist, the possible values of parameters ν are restricted.
The proof of the corollary is a special case, where one does not have to introduce the integrals over t. The existence of the inverse A −1 is ensured by theorem 1.4.
Remark: A more general version of this theorem can be found in [3, thm. 5.6] with an independent proof. In [72] a similar technique is used to obtain Mellin-Barnes representations from Feynman integrals.
The convergence of those Mellin-Barnes representation is discussed in [34] and general aspects of multivariate Mellin-Barnes integrals can be found in [58].
The corollary can alternatively be proven by splitting the term G −ν 0 by the multinomial theorem in a multidimensional power series and solve the Mellin transform of this power series by a generalized version of Ramanujan's master theorem [32].
The type of the integral which appears in equation (1.10) is also known as multivariate Fox H-function [34] and the connection between Feynman integrals and Fox H-function was studied before [15,39,40].
Since the number of monomials in the Symanzik polynomials increases fast for more complex Feynman graphs, the Mellin-Barnes representation of theorem 1.2 does not provide an efficient way to calculate those integrals. One exception are the massless 2-point functions consisting in n lines and having the loop number L = n − 1. These so-called "banana graphs" or "sunset-like" graphs, have only n + 1 monomials in G and satisfy therefore the condition in the corollary.
In order to illustrate the theorem 1.2 about Mellin-Barnes representations, we finish this section with a simple example.  Choosing σ = {1, 2, 3}, the corresponding Feynman integral in the Mellin-Barnes representation of theorem 1.2 is given by (1.16) For the correct contour prescription the poles have to be separated such that there exist values δ satisfying max{0, −ν 0 + ν 1 } < δ < min{ν 0 − ν 2 , −ν 0 + ν 1 + ν 2 }. In the formulation to be introduced in the following section 1.3, this is equivalent to claim a full dimensional Newton polytope of G. In this case, by Cauchy's theorem the integral evaluates simply to a Gaussian hypergeometric function After restoring the original prefactors and coefficients z 1 = z 2 = 1, z 3 = p 2 − m 2 1 , z 4 = m 2 1 and ν 0 = d 2 it agrees with the expected result (1.18)

Mathematical Interlude: Convex Polytopes
In the following, it will be fruitful to describe properties and solutions of Feynman integrals from a perspective of convex polytopes. In order to clarify some terminology, we give a short review of the basic concepts of convex polytopes. Readers which are familiar with convex polytopes can skip this interlude. More details about polytopes, customized to the hypergeometric approach, can be found in the appendix. For further treatments we refer to [11,13,21,35].
A convex polytope P is defined as the convex hull of a finite set of points A = {a 1 , . . . , a N } where a i ∈ R n . Additionally, if the points a i ∈ Z n form an integer lattice, P is called a convex lattice polytope. Since in this discussion all polytopes will be convex polytopes, we will call them simply "polytopes" in the following. As a fundamental result of polytope theory, every polytope can also be written as a bounded intersection of half-spaces where b ∈ R k and M ∈ R k×n are real, m T j denotes the rows of the matrix M and · is the standard scalar product.
A subset of P having the form F := {µ ∈ P |c µ = β}, with c ∈ R k×n and where every point of the polytope µ ∈ P satisfies the inequality c µ ≤ β, is called a face of P . Faces of dimension zero are called vertices and faces of dimension n − 1 are called facets. Thus, if equation (1.20) consist in a minimal set of inequalities, the facets of P are given by m T j µ = b j . The polytope without its faces is called the relative interior relint P of the polytope P .
The dimension of a polytope P = Conv(A) is defined as the dimension of its affine hull, and thus one can easily see As before in equation (1.8), we will denote the matrix A with an additional row (1, . . . , 1) as A = 1 A . If a polytope P ⊂ R n has the dimension n it is called to be full dimensional and degenerated otherwise. For full dimensional polytopes it is meaningful to introduce a volume of polytopes. As the volume vol 0 (P ) ∈ N of a lattice polytope we understand a volume, which is normalized such that the standard simplex has volume 1. Therefore, this volume is connected to the standard Euclidean volume vol(P ) by a faculty of the dimension vol 0 (P ) = n! vol(P ). For a full dimensional simplex P = Conv(A) the volume is given by the determinant vol 0 (P ) = | det A|. Thus, one way of calculating the volume of a polytope is by dividing the polytope in simplices.
A subdivision of a polytope P into simplices, where the union of all simplices is the full polytope and the intersection of two distinct simplices is a proper face of the polytope (including the empty set), is called a triangulation.
A triangulation T (ω) = {σ 1 , . . . , σ r } of a polytope P = Conv(A) is called regular, if there exists an height vector ω ∈ R N , such that for every simplex σ i of this triangulation there exists another vector r i ∈ R n+1 satisfying where · denotes the scalar product. It can be shown, that every convex polytope admits always a regular triangulation [21]. If all σ ∈ T (ω) belongs to simplices with volume 1 (i.e. | det A σ | = 1) the triangulation is called unimodular.
In order to reduce later more complicated Feynman integrals to simpler one, we introduce subtriangulations.

Properties of Feynman Integrals
In the definitions of the Feynman integrals we have omitted to discuss the convergence of those integrals. E.g. by power counting [76] the convergence behavior of the original Feynman integrals in momentum space is well known, as well as for the parametric Feynman integral (1.3) [57]. Therefore, we discuss the convergence shortly for the representation (1.5) and (1.9), respectively. This discussion involves the perspective of polytopes, which allows a clear and short notation. The theorems are mostly direct implications of the work of [3,54,65] and proofs can be found there.
In the description of (1.20) this is equivalent to demanding b j Re ν 0 − m T j · Re ν > 0 for 1 ≤ j ≤ k. Furthermore, if the Newton polytope ∆ G is not full-dimensional, the Feynman integral does not converge absolutely for any choice of ν 0 and ν.
The second statement means, that Feynman integrals which do not have a fulldimensional Newton polytope ∆ G are neither dimensionally nor analytically regularizable. This result is not surprising, since if the Newton polytope is not full dimensional, the polynomial G has a special homogeneous property. A polynomial G is homogeneous in such that way, that there exists numbers c 0 , . . . , c n ∈ Z not all zero, such that holds. If G corresponds to a massless tadpole graph, G is homogeneous in this sense, since massless tadpole graphs have F = 0 and the first Symanzik polynomial is homogeneous of degree L. Thus, there will be no values for ν = (ν 0 , ν) ∈ C n+1 that those integrals converge. However one usually sets those tadpole graphs to zero by including a counterterm in the renormalization procedure, which removes all tadpole graphs [70] (see also [6, sec.

3.1]).
Since the Feynman integral in (1.9) is only equivalent to the original Feynman integral in the meromorphic extension, the convergence region of (1.9) has no deeper physical meaning. Using the following theorem one can deduce the meromorphic continuation of (1.9) to the whole complex plane, which determines the physical relevant poles of the Feynman integral.
Then one can rewrite the Feynman integral as where Φ A (ν, z) is an entire function with respect to ν ∈ C n+1 . Example 1.3: Consider the example from above which correspond to figure 1. For the relative interior of the Newton polytope one obtains from the facet representation, the region of convergence (with Re ν 0 > 0) which enables us to separate the poles of the Feynman integral in the Γ functions This result is remarkable in many different ways. Firstly, it guarantees that the Feynman integral can be meromorphically continued to the whole complex plane, which confirms dimensional and analytical regularization. Secondly, equation (1.25) gives an easy method to calculate the (possible) poles of the Feynman integral. And thirdly, in the expansion one can focus only on a Taylor expansion of Φ A instead of a Laurent expansion of J A . Thus, one can determine the coefficients by differentiating, which makes the procedure much easier.
These theorems do not rely on any special properties of the Symanzik polynomials. In fact, Feynman integrals are just a subset of Euler-Mellin integrals. The following lemma is a simple implication from the properties of Symanzik polynomials to be at most quadratic, and to be homogeneous of degree L and L + 1, respectively. c) The Newton polytope ∆ G , corresponding to a Feynman graph, has no interior points.
Thus, we have reformulated the Feynman integral as an Euler-Mellin integral, which defines a meromorphic function in ν ∈ C n+1 . To consider the Euclidean region we have restricted 3 the discussion to the right half space z ∈ C N with Re z j > 0. Further, in corollary 1.3 we had a class of Feynman integrals which provides a simple and analytic solution. These integrals will be helpful to find boundary values for the following partial differential equation systems. The Newton polytopes ∆ G , which define the convergence regions and determine the whole structure of the Feynman graphs, which is necessary to evaluate Feynman amplitudes, are relatively well behaved.

General Hypergeometric Functions
Since the first hypergeometric function was studied by Euler and Gauss more than 200 years ago, many different generalization of hypergeometric functions were introduced: Pochhammer series p F q , Appell's, Lauricella's and Kampé-de-Fériet functions, to name a few. Those functions can be characterized in three different ways: by series representations, by integral representations and as solutions of partial differential equations. Therefore, there are in principle three different branches to generalize the notion of a hypergeometric function.
The most general series representation goes back to Horn [36] and was later investigated by Ore and Sato (a summarizing discussion can be found in [26]). A Horn hypergeometric series is a multivariate power series in the variables where ratios of coefficients c(k+e i ) are rational functions in k 1 , . . . k r and where e i is the standard basis in the Euclidean space. Thus, the coefficients can be represented mainly by a product of Pochhammer symbols 4 , which are defined as (a) n := Γ(a + n) Γ(a) .
Negative integers in arguments of Γ functions can be avoided by considering appropriate limits. In the terms of the series of Horn hypergeometric functions we consider a ∈ C as complex numbers and n as integer combinations of the k 1 , . . . , k r . Among many beautiful properties, derivatives with respect to the parameters of Horn hypergeometric are again Horn hypergeometric functions [16]. For further studies of Horn hypergeometric functions we refer e.g. to [61].

GKZ Hypergeometric Functions
Since the late 1980s the theory of general hypergeometric functions was reinvented by the characterization of the partial differential equation system by Gelfand, Graev, Kapranov, Zelevinsky and collaborators [26][27][28][29][30][31]. This approach combines the different characterizations of hypergeometric functions and gives a comprehensive method to analyze and describe general hypergeometric functions. In this section we introduce the Gelfand-Kapranov-Zelevinsky (GKZ) hypergeometric system, which is a system of partial differential equations and discuss roughly how to solve this system by power series. The following recapitulation is adapted for the application to Feynman integrals and will miss some generality in order to keep the discussion short. Since the theory of general hypergeometric functions involves many mathematical aspects inter alia algebraic geometry, combinatorics, number theory and Hodge theory, we refer for more detailed studies to [2,19,26,31,62,71,75].
Let A ∈ Z (n+1)×N be an integer matrix with n+1 ≤ N , rank A = n+1 and including a row of the form (1, 1, . . . , 1). The ladder condition means that A lies in an n-dimensional affine hyperplane of Z n+1 . Without loss of generality we can consider that A is of the form A = 1 A as in equation (1.8). Then the Gelfand-Kapranov-Zelevinsky (GKZ) hypergeometric system is defined as the D-module, consisting of toric and homogeneous differential operators is the Euler operator and β ∈ C n+1 . By D we denote the ideal generated by components of D over C. Solutions Φ(z) of these differential systems One of the most important properties of those hypergeometric systems (2.3) is to be holonomic, i.e. the dimension of the solution space is finite and we can give an appropriate basis of the solution space. There are several possibilities to construct those bases. In the following we discuss a solution in terms of multivariate power series, which was the first solution invented in [28]. Furthermore, there are solutions in terms of different sorts of integrals [49]. The following discussion is simplified for the later approach and based on [23] and [71].
Let A ∈ Z (n+1)×N be an integer matrix with n + 1 ≤ N and full rank as before. Furthermore consider the corresponding lattice L := ker Z A = {(l 1 , . . . , l N ) ∈ Z N |l 1 a 1 + . . . + l N a N = 0}, which has rank L = N − n − 1 =: r by the rank-nullity theorem. Then for ξ ∈ C N the formal series is called Γ -series. It turns out, that these Γ-series are formal solutions of the GKZ system (2.3) for an appropriate choice of ξ: Lemma 2.1 [Γ-series as formal solutions of GKZ hypergeometric systems [28,71]]: Let L be the corresponding lattice to A and ξ ∈ C N satisfying Aξ + β = 0. Then the series z r−u (with an appropriate limit, respectively). Furthermore one can add an element of L to ξ, without changing the which shows that the Γ-series satisfies the toric equations. For the homogeneous equations one considers The restriction Aξ + β = 0 allows in general many choices of ξ. Let σ ⊆ {1, . . . , N } be an index set with cardinality n + 1, such that the matrix A restricted to columns of that index set σ is invertible, det A σ = 0. Due to the assumption rank A = n + 1 those index sets always exist. Denote byσ = {1, . . . , N } \ σ the complement of σ. If one sets ξ σ = −A −1 σ (β + Aσk) and ξσ = k the condition Aξ + β = 0 is satisfied for any k ∈ C r . On the other hand we can split the lattice L = {l ∈ Z N |Al = 0} in the same way A σ l σ + Aσlσ = 0 and obtain a series only over lσ In order to simplify the series one can choose 5 k ∈ N r 0 , since terms with (k + lσ) i ∈ Z <0 will vanish. The Γ-series depends now on k and σ Also the choice k ∈ Z r would be possible, but it does not change the series, see [23, lemma 3.2.] where Λ k = {k + lσ ∈ N r 0 |Aσlσ = ZA σ } ⊆ N r 0 for any k ∈ N r 0 . Therefore, the Γ-series is turned into a power series.
Remark: In the unimodular case | det A σ | = 1, the coefficient matrix A −1 σ Aσ ∈ Z (n+1)×r is an integer matrix and therefore it is . In order to show that Γ-series are actual solutions of the GKZ system and not only formal ones, one has to prove that Γ-series converge for some z ∈ C N . By an application of the Stirling approximation it can be shown, that the Γ-series always converge absolutely for sufficiently small values of the variables x j := . A proof of the absolute convergence of Γ-series can be found in the appendix 6.2.
Another issue is also that the Γ-series can be identical to zero, which is also inconvenient in order to construct a solution space. The Γ-series is zero for all z ∈ C N , if and only if for all λ ∈ Λ k the expression A −1 σ (β + Aσλ) contains at least one positive integer entry. To avoid these cases one considers generic β ∈ C n+1 : is a strictly positive integer for all λ ∈ N r 0 one says that β is very generic with respect to σ. In the unimodular case this is equivalent to claim, that for the components i which satisfy ( Thus, typically non generic cases arise for even integer dimensions, which should not be surprising, since for these dimensions the Feynman integrals mostly diverge, which is the reason why we exclude this values already.
To normalize the first term of the power series to 1, we will deal in the following with a slightly different version of the Γ-series Aσλ denotes the (multivariate) Pochhammer symbol (a) n := j Γ(a j +n j ) Γ(a j ) . Remark: In the unimodular case | det A σ | = 1 one can rewrite the Γ-series by Pochhammer identities.
As mentioned above the holonomic rank is finite. For very generic β one can determine the holonomic rank by a polytope corresponding to A = 1 A .
Theorem 2.2 [Holonomic rank of GKZ systems [19,29,30]]: Consider a GKZ system M A (β) with arbitrary A ∈ Z (n+1)×N and very generic β ∈ C n+1 . Let Conv(A) be the corresponding convex polytope and denote by vol 0 the normalized Euclidean volume, such that the standard simplex has a volume equal to 1. Then the holonomic rank of the GKZ system is equal to the volume of the polytope Conv(A) rank M A (β) = vol 0 (Conv(A)) . (2.11) This means that one needs vol 0 (Conv(A)) linearly independent solutions to construct the solution space. The regular triangulations of the polytope Conv(A) provide a construction of linearly independent Γ-series. Definitions and basic properties of the regular triangulations can be found in the mathematical interlude of section 1.2 and in the appendix.
In the following we will only discuss the case of unimodular triangulation, since almost all Feynman integrals admit an unimodular triangulation. We will motivate this restriction in section 3.5 in more detail. Nevertheless, there is a simple generalization to the non-unimodular case, which can be found e.g. in [23], [49]. [23,29]]: Let T be a regular unimodular triangulation and let β be very generic with respect to every σ ∈ T . Then the Γ-series {ϕ σ } σ∈T form a basis of the solution space of the hypergeometric GKZ system M A (β). Furthermore, all these Γ-series have a common region of convergence.

Theorem 2.3 [Solution Space of GKZ
To conclude, we defined holonomic systems of partial differential equations, which can be characterized by a matrix A ∈ Z (n+1)×N and a vector β ∈ C n+1 . Furthermore, for generic values of β we are able to construct the whole solution space in terms of power series by regular triangulations of the polytope Conv(A).

Feynman Integrals as Hypergeometric Functions
It is one of the first observations in the calculation of simple Feynman amplitudes, that Feynman integrals evaluate mostly to hypergeometric functions. This observation was leading Regge to the conjecture that Feynman integrals are always hypergeometric functions and he based his conjecture on the partial differential equations which are satisfied by the Feynman integral [60].
Typically, those hypergeometric functions also appear in the often used Mellin-Barnes approach. This is a consequence of Mellin-Barnes representations with integrands consisting in a product of Γ functions, which can be identified by the hypergeometric Fox H-functions [15,39,40] or which can be evaluated to some series-based hypergeometric functions, like the Appell or Lauricella functions by application of the residue theorem [8,24]. But except of some special cases, like the one-loop integrals [24], this correspondence is more or less unproved, which is also due to the fact that multivariate Mellin-Barnes integrals can be highly non trivial [58].
A new opportunity to examine the correspondence between hypergeometric functions and Feynman integrals is the Gelfand-Kapranov-Zelevinsky approach. It was already stated by Gelfand himself, that "practically all integrals which arise in quantum field theory" [30] can be treated with this approach. Recently, the connection between Feynman integrals and A-hypergeometric function was concertized in [52] and, independently of the present paper, in [20].
Based on the Lee-Pomeransky representation of Feynman integrals it is a standard application to show that generic Feynman integrals are A-hypergeometric, analogue to the examples in [71]. Thus the generic Feynman integral is an A-hypergeometric function.
Proof. Firstly, we show that the generic Feynman integrals satisfy the toric part {∂ u − ∂ v |Au = Av, u, v ∈ N N }. Derivatives of the Feynman integral with respect to z result in where |u| := i u i . From the row (1, 1, . . . , 1) in A it follows immediately that |u| = |v|. Therefore, one obtains the same equation for v. Secondly, consider for the homogeneous part Aθ A derivative with respect to s completes the proof for s = 1.
Thus, as suggested already by Gelfand and confirmed in [52] and [20], every generic Feynman integral with an Euclidean region satisfies the GKZ hypergeometric system and can be treated within the framework of GKZ. This will allow inter alia a series representation of the Feynman integral.

Hypergeometric Series Representations of Generic Feynman Integrals
As stated in theorem 3.1 generic Feynman integrals are A-hypergeometric functions. Thus, one can directly apply the results of hypergeometric GKZ systems from section 2. Consider first the case ν ∈ (C \ Z) n+1 in order to satisfy the very genericity of ν and to be sure that the Feynman integral has no poles. Later we can relax this strict claim. We fix a regular, unimodular triangulation T of the Newton polytope ∆ G of the Lee-Pomeransky polynomial G = U + F . Then by theorem 2.3, we can write the generic Feynman integral as a linear combination of Γ-series where the z ∈ C N are defined in the positive half-space Re z i ≥ 0 in the region where the Γseries converge. Thus, one has to determine the meromorphic functions C σ (ν) in order to get a series representation of Feynman integrals. This can be done by comparing equation (3.4) with boundary values of the Feynman integral. As we will see in the following, the Feynman integrals transmit their functions C σ (ν) to simpler Feynman integrals. This enables us to reduce Feynman integrals to the case described in corollary 1.3 and derive an analytic expression of the functions C σ (ν). For this purpose, consider a generic Feynman integral J A with a Newton polytope P = ∆ G = Conv(A), consisting of the vertices a 1 , . . . , a N , and an unimodular triangulation T = {σ 1 , . . . , σ r , η 1 , . . . , η s } of P . Furthermore, let T = {σ 1 , . . . , σ r } be a proper subtriangulation of T . Denote the vertices of the convex polytope P , which correspond to the subtriangulation T , by a 1 , . . . , a M with M < N and the corresponding Feynman integral 7 by J A . Thus, in the second Feynman integral some monomials in G z are missing Applying the results from the previous sections to both Feynman integrals independently one obtains on the one hand with the Γ-series and on the other hand Due to the construction it is a M +1 , . . . , a N ∈ r i=1σ i . Therefore, in the limit z M +1 , . . . , z N → 0 only some zσ will be affected and it is simply Apart from that, as a consequence of the homogeneous differential equations Aθ + ν all solutions of a GKZ hypergeometric system have to satisfy the scaling property φ(ν, sz) = s −ν 0 φ(ν, z). It can be easily seen that Γ-series satisfy this property. In the power series part of a Γ-series the variables z only appear as ratios, such that they are scaling invariant (in accordance with lemma 6.3). However, the monomial z of the power series will give the scaling property since s i∈η (−A −1 η ν) i = s −ν 0 according to lemma 6.3.
Since in every Γ-series φ η i some of the variables z M +1 . . . , z N are contained in the monomial, the scaling property will be violated in the limit z M +1 . . . , z N → 0. Thus, the functions lim z M +1 ,...z N →0 φ η i can not be linearly dependent of {φ σ i }. Since one already has r linearly independent solutions φ σ i and therewith has a full-dimensional solution space for the Feynman integral J A , the Γ-series φ η i have to vanish lim z M +1 ,...,z N →0 φ η i (ν, z 1 , . . . , z N ) = 0 . (3.12) Thus, the GKZ systems behave naturally, as also mentioned in [26]: If one deletes a vertex of Conv(A), the Γ-series which correspond to simplices containing this vertex will vanish. This leads to a very simple connection between subtriangulations of Feynman integrals. Applying the limit z M +1 , . . . , z N → 0, the meromorphic functions C σ (ν) will not be affected and one obtains C σ i (ν) = C σ i (ν).
Thus, one can determine the meromorphic functions C σ (ν) by considering simpler Feynman integrals which refer to subtriangulations, where by a simpler Feynman integral we mean a Feynman integral where the Lee-Pomeransky polynomial G has less monomials. In this way one can define ancestors and descendants of Feynman integrals by deleting or adding monomials to the Lee-Pomeransky polynomial G. E.g. the massless one-loop bubble graph is a descendant of the one-loop bubble graph with one mass, which itself is a descendant of the full massive one-loop bubble. Those ancestors and descendants do not necessarily correspond to Feynman integrals in the original sense, since one can also consider polynomials G which are not connected to graph polynomials anymore.
As a trivial subtriangulation of an arbitrary triangulation one can choose one of its simplices. In doing so, one can relate the prefactors C σ (ν) to the problem where only one simplex is involved. For such problems, one can solve the Feynman integral easily as seen in corollary 1.3. Therefore, for an unimodular triangulation one can find that the prefactors are simply given by which results in the following theorem:

Theorem 3.2 [Series representation of Feynman integrals]:
Let T be a regular, unimodular triangulation of the Newton polytope ∆ G = Conv(A) corresponding to a generic Feynman integral J A . Then the generic Feynman integral can be written as where the series have a common region of convergence. This representation holds for generic ν ∈ C n+1 , which means that ν has to be chosen such that the Feynman integral has no poles and none of the power series in (3.14) will be identical to zero.
In the physical relevant limit z → (1, 1, m 2 1 − p 2 , m 2 1 ) and ν = (2 − , 1, 1) one can easily evaluate the series which agrees with the expected result. The series representation which can be obtained by the triangulation T 2 , as well as the former result in example 1.1, are equivalent to this result by transformation rules of the 2 F 1 function.
Thus, we found a series representation for generic Feynman integrals which admit an unimodular triangulation. Based on experience, it is reasonable to conjecture that every off-shell Feynman graph admits an unimodular triangulation. But also for nonunimodular triangulations, the Feynman integral can be referred back to a case which admits an unimodular triangulation. We will examine the case of non-unimodular triangulations and this relation in section 3.5. Typically, a Feynman graph admits many different possibilities to triangulate its corresponding Newton polytope. Therefore, one usually obtains a large number of series representations. This is not surprising, since hypergeometric functions satisfy many transformation formulas and can be converted to other hypergeometric functions.
Therefore, in practical computations one can choose a series representation, which converges fast for the given kinematics and evaluate the Feynman integral numerically by considering the first summands of every series.
These series representations -like the whole GKZ approach -are only valid for generic Feynman integrals. In principle this series representation is also true for fixed values of z. But for non-generic values of z it can happen, that the series in (3.14) do not converge anymore. We address this point in the following.

The Limit from Generic Feynman Integrals to Non-Generic Feynman Integrals
Up to this point, the theorems were statements about the generic Feynman integral. The next natural question is how one can perform the limit to the non-generic, "ordinary" Feynman integral, i.e. the limit, where the extra variables, which are connected to the first Symanzik polynomial are set equal to 1. In this limit the convergence behavior of the Γ-series can be changed. Consider a region D ⊆ C n+1 where the Feynman integral (1.5) converges for ν ∈ D. Since in this case the linear combination (3.4) is still finite, there can arise only two problems: a) every series converges separately, but they do not have a common convergence region anymore or b) some of the Γ-series diverge, but the linear combination is still finite. In the first case a) the convergence criteria for the variables of the Γ-series x j = (zσ) j i (z σ ) −(A −1 σ Aσ) ij i exclude each other for different σ ∈ T . In the second case b) the variables x j become constants after the limit (usually equals 1), which can be outside of the convergence region.
Because of these possible issues, it can be difficult to perform the limit. Fortunately there is a strategy to tackle these problems in many cases by transformation formulas of hypergeometric functions. E.g. for the 2 F 1 hypergeometric function, there is a well known transformation formula [56] which can be applied to change a limit x j → 1 to the much simpler case of a limit x j → 0. We illustrate this method with an example. For the 2-loop sunset graph with two different masses, inter alia there appears the hypergeometric series where one has to consider the limit (z 1 , z 2 , z 3 , z 4 , z 5 , z 6 , z 7 , z 8 ) → (1, 1, 1, m 2 2 , m 2 1 + m 2 2 − p 2 1 , m 2 1 , m 2 2 , m 2 1 ). In this limit it appears the term (−1) k 4 , which is not in the convergence region for small values of > 0 anymore. Therefore, we evaluate the k 4 series carefully and write With the transformation formula for the 2 F 1 function, one can split the series in a convergent and a divergent part Comparing the divergent part with the other Γ-series, which occurs in the calculation of the sunset graph with two masses, one can find another divergent series which exactly cancels this divergence. This cancellation has always to happen, since the linear combination has to be finite. Therefore, one can derive a convergent series representation also for non-generic Feynman integrals by considering convenient transformation formulas of hypergeometric functions. In principle, this procedure will work in general, but it can be necessary to involve more complicated transformation formulas than the well-studied 2 F 1 transformation and also the cancellation of the divergences may be not so obvious as in this case. Fortunately, many graphs (also for L > 1) can even be evaluated with the transformation formula of the 2 F 1 function.
In fact, this limit can reduce the dimension of the solution space, which is the expected behavior. For generic variables z ∈ C N the dimension of the solution space is equal to vol 0 ∆ G according to theorem 2.2. In contrast for non-generic values of z ∈ C N the dimension of the solution space is equal to the Euler characteristic (−1) n χ((C ) n \ {G = 0}) [6], which is in general smaller or equal to the volume of the Newton polytope (−1) n χ((C ) n \ {G = 0}) ≤ vol 0 ∆ G . (3.22) Thus, by calculating the Euler characteristic we can count the expected dependencies in the limit from generic to non-generic Feynman integrals. According to [6] it is meaningful to define the number of master integrals as the dimension of the solution space. Therefore, one can see the linear combination of Γ-series (3.4) as a procedure similar to a decomposition of a general Feynman integral into a basis of master integrals. This analogy mirrors also the existence of different decompositions (which corresponds to different triangulations) and their transformation into each other by shift relations [47].

The Expansion of Hypergeometric Series Representations
In the calculus of dimensional regularization [73] one is usually interested in the Laurent expansion of the Feynman integrals around = 0 where d = 4−2 . Due to the theorem 1.5 one can relate this task to the Taylor expansion of the hypergeometric series representation. Thus, one has simply to differentiate the Horn hypergeometric series. As pointed out in [16], the derivatives with respect to parameters of Horn hypergeometric series are again Horn hypergeometric series of higher degree. By the identities (a) m+n = (a) m (a + m) n and (a) rn = r rn r−1 j=0 a+j r n for r ∈ Z >0 one can reduce all derivatives to two cases [16] Thus, Horn hypergeometric functions do not only appear as solutions of Feynman integrals with unimodular triangulations, but also in every coefficient of the Laurent expansion of those Feynman integrals. Therefore, the class of Horn hypergeometric functions is sufficient to describe almost all Feynman integrals and their Laurent expansion. Hence, it is not surprising that also the combinatorial structure of Feynman integrals is reflected in the Horn hypergeometric functions. For instance, the relations between different Feynman integrals can be derived by transformation formulas of hypergeometric functions and vice versa [47]. Another way to expand the Horn hypergeometric functions around = 0 in some cases could be the approach of S-and Z-sums [50], which are related to multiple polylogarithms and related functions. Unfortunately, many examples, which can be generated by the GKZ approach, belong not to the known algorithms given in [50].

Applications to other Parametric Representations
The GKZ mechanism with the resulting series representation used above is not limited to the Feynman integral in the Lee-Pomeransky representation. It is a method which can be used for all integrals of Euler-Mellin type including one polynomial. Thus, there are more applications in the Feynman integral calculus. For example the Feynman parametric representation (1.3) is of that form, if either the first Symanzik polynomial U or the second Symanzik polynomial F drops out. This is the case e.g. for all 1-loop Feynman integrals, where the first Symanzik polynomial U evaluates to 1 after applying the δ-function. With ν n = 1 after one integration the Feynman integral is of the form of an Euler-Mellin integral we denotex = (x 1 , . . . , x n−1 ),ν = (ν 1 , . . . , ν n−1 ) andF = F | xn→1− n−1 i=1 x i . Also for so called "marginal" Feynman integrals [9], where ω = d/2, the first Symanzik polynomial drops out of the representation (1.3) and one obtains the integral In contrast the second Symanzik polynomial drops out for ω = 0, which is the wellstudied case of periods (e.g. [12]) Note that this integral is highly non-generic from the perspective of the GKZ approach. Thus, the case where the second Symanzik polynomial F remains is usually much easier, since one has not to introduce extra variables in the first Symanzik polynomial U .
Last but not least, also the Baikov representation [6,33] is a possible candidate to apply the GKZ approach as well.

Non-Unimodular Triangulations of Feynman Polytopes
The treatments above were specialized to unimodular triangulations only. This strategy has various reasons. Without much effort, one could extend theorem 2.3 also to the case of non-unimodular triangulations and write the Feynman integral as a linear combination of Γ-series as in equation (3.4). In contrast, for a non-unimodular triangulation one can not determine the meromorphic functions C σ (ν) such as easy as in theorem 3.2. Nevertheless, one can reduce Feynman integrals without unimodular triangulations to subtriangulations as described in section 3.1.
However, after checking common Feynman integrals up to three loops, we are leaded to the conjecture that all off-shell Feynman graphs admit at least one unimodular triangu-lation. This conjecture seems also likely by considering the very specific form of Newton polytopes ∆ G arising in Feynman integrals (lemma 6.2).
In contrast some on-shell graphs do not admit an unimodular triangulation. For instance the on-shell full massive 1-loop bubble with G = x 1 + x 2 + m 2 1 x 2 1 + m 2 2 x 2 2 does not allow an unimodular triangulation. However, since the off-shell fully massive 1-loop bubble admits unimodular triangulations, one can treat the on-shell Feynman integral as a limit of the off-shell version.
This behavior holds in general: one can always add monomials to the Lee-Pomeransky polynomial G → G , such that the Newton polytope ∆ G admits an unimodular triangulation. This can be seen by a result of Knudsen et al. [46]: 13,46]]: For every lattice polytope P there is an integer k ∈ N such that the dilated polytope kP := {kµ|µ ∈ P } admits an unimodular triangulation.
Thus, if the Newton polytope ∆ G of a Feynman integral does not admit an unimodular triangulation, one can dilate the polytope ∆ G until the polytope k∆ G admits an unimodular triangulation. Then the original Feynman integral can be obtained as a limit where all additional vertices vanish and where we scale the propagator powers by k and the whole integral by k n Thus, one can always write the Feynman integral as a limit of another integral where its Newton polytope allows an unimodular triangulation. That is the reason, why we can focus only on the unimodular triangulations and vindicates the procedure above.

Advanced Example: Full Massive Sunset
To illustrate the GKZ method stated above, as well as to show the power of this approach, we calculate a series representation of the sunset Feynman integral with three different masses according to figure 3. The corresponding Feynman graph consists in n = 3 edges and the Lee-Pomeransky polynomial includes N = 10 monomials The rank of the kernel of A is equal to r = N − n − 1 = 6 and therefore we will expect 6-dimensional Γ-series. Moreover, the polytope ∆ G = Conv(A) has the volume vol 0 (Conv(A)) = 10 (calculated with polymake [25]), which leads to 10  in order to get series, which converge fast for highly relativistic kinematics m 2 i m 2 1 + m 2 2 + m 2 3 − p 2 . Further, we set ν i = 1 and d = 4 − 2 . In the limit z → (1, 1, 1, m 2 3 , m 2 3 , m 2 2 , m 2 1 + m 2 2 + m 2 3 − p 2 1 , m 2 1 , m 2 2 , m 2 1 ) the series φ 1 , φ 3 , φ 5 , φ 6 , φ 8 and φ 9 are divergent for small values of > 0. By the method described in section 3.2 one can split all these series by the transformation formula for the 2 F 1 hypergeometric function in a convergent and a divergent part. The divergent parts of these series cancel each other. In doing so the resulting Γ-series have linear dependences and the dimension of the solution space will reduce from 10 to 7.
By applying all these steps one arrives at the following series representation of the full massive sunset integral where the Γ-series are given by 3 −p 2 and s = m 2 1 +m 2 2 +m 2 3 −p 2 . All these series converge for small values of x i and the series representation can be obtained by a very simple algorithm, which is a straightforward implementation of the steps described in section 3. In fact, some of these Γ-series are related to each other. One can reduce the whole system only to φ 1 , φ 2 , φ 5 and φ 10 by the relations φ 1 (x 1 , x 2 , . By these relations one can also verify the expected symmetry of the Feynman integral under the In order to expand the Feynman integral I A for small values of > 0 one can use the methods described in section 3.3 or alternatively by expanding each Γ-function separately. The latter requires to distinguish between positive and negative integers in the argument of the Γ-function (4.7) By the distinction of cases between positive and negative arguments in the Γ functions many terms arise, which are easily manageable by a CAS but which are space-consuming in print, wherefore we omit to state these results here. The correctness of these results was checked numerically by FIESTA [67] with arbitrary kinematics and masses, satisfying x i < 0.5. For small values of x i the resulting series converges fast, such that for a good approximation one only has to take the first summands into account. An upper bound for the errors in a finite summation can be estimated by a majorant geometric series, similar to the procedure in theorem 6.6. Vice versa one can determine the number of required summands for a given error bound. Furthermore, the summands of a Horn hypergeometric series always have a rational ratio. Thus, in a numerical calculation one only has to evaluate rational functions. In this way Horn hypergeometric series have a relatively simple and controllable numerical behavior in the case |x i | 1.

Conclusion and Outlook
We showed in this article that Feynman integrals can be described as hypergeometric functions. Namely we showed that a) every generic Feynman integral is an A-hypergeometric function, that b) every generic Feynman integral which admits a regular, unimodular triangulation has a representation in Horn hypergeometric functions and that c) all scalar, dimensional regularizable and Euclidean Feynman integrals can be written at least as a limit of a linear combination of Horn hypergeometric functions. Furthermore, the latter (c) is also true for all coefficients in a Laurent expansion of the Feynman integral in a dimensional or analytic regularization. Since hypergeometric functions are mostly represented in terms of integrals including polynomials, Mellin-Barnes integrals including Γ-functions or series including Pochhammer symbols, it is not surprising that also Feynman integrals appear normally in one of those representations. From the perspective of general hypergeometric functions the common ground of these representations is an integer matrix A ∈ Z (n+1)×N or equivalently a Newton polytope ∆ G = Conv(A). By the triangulation of this Newton polytope we derived an analytic formula for a hypergeometric series representation of the Feynman integral. As there are in general many different ways to triangulate a polytope, there are also many different series representations possible. Similar to the Feynman integral those hypergeometric series satisfy different relations between each other and can therefore also be transformed in equivalent representations. However, the series representations can differ in their convergence behavior. Thus, there can be series representations which converge fast for given kinematics so that we sometimes only need the first summands in order to reach a high precision numerical approximation.
For the purpose of a practical usage of this concept, we discussed possible obstacles which can appear in the concrete evaluation and gave some strategies to solve them. The procedure described in section 3 and illustrated in section 4 works similarly in all cases and can be read as an algorithm.
Besides numerical applications, there are structurally interesting implications for the Feynman integral. Since in both subjects similar questions appear, the hypergeometric perspective opens new ways to analyze the Feynman integral. For instance the singularities of hypergeometric functions, and therefore simultaneously the singularities of Feynman integrals, are given by the A-discriminants and the coamoebas.
The series representation given in the present work, is only one of many consequences that we can obtain by the connection between Feynman integrals and hypergeometric functions. On the one hand we outlined in section 3.4 that we also can apply the GKZ approach to other integrals appearing in the Feynman calculus. On the other hand one can also solve the GKZ system with other functions to obtain representations of the Feynman integral e.g. in terms of other integrals. This leads us to the question, if one can relate Feynman integrals with the GKZ approach to other Feynman integrals in a consistent way, in order to get linear relations between Feynman integrals similar to the IBP approach.
Furthermore, it would be desirable to get a better understanding of the relations between Feynman integrals and their diagrams sharing parts of the matrices A. The simplest connection between those Feynman integrals are the ancestors and descendants of Feynman integrals mentioned in section 3.1. Lastly, we leave several technical questions open to solve in future work. For instance, a proof of the conjecture that off-shell Feynman integrals always admit unimodular triangulations together with a better understanding of the behavior in non-generic cases as described in section 3.2, would show that all off-shell Feynman integrals can be written as Horn hypergeometric functions and not only as a limit of Horn hypergeometric functions.
In this article we suggested a new perspective of an old idea to characterize Feynman integrals. This hypergeometric perspective can be useful for practical computation but also for a better conceptual understanding of Feynman integrals. Still, this is merely a starting point for future research on the correspondence between Feynman integrals and hypergeometric functions in the GKZ theory.

Convex Polytopes & Triangulations
Convex polytopes constitute one of the main concepts in the present paper and are the key link between Feynman integrals and hypergeometric functions. They characterize the whole structure of a Feynman graph, which is necessary to evaluate the Feynman integral as a function of kinematics and masses. Furthermore they determine the region of convergence of the Feynman integral (theorem 1.4), the poles and the meromorphic continuation of the Feynman integral (theorem 1.5) and their triangulations give series representations (theorem 3.2) in terms of Horn hypergeometric series.
In section 1.2 we already invented the most important concepts which are essential for this approach. At this point, we give some more details which are closely related to this approach and which give an illustrative viewpoint on polytopes. For further treatments we refer to [11,13,21,35].
In the mathematical interlude in section 1.2 we introduced the two ways to represent a polytope (by vertices and by the intersection of half-spaces). For general polytopes the transformation of these two representations can be complicated. However for simplices there is a simple connection. Proof. Since P is a full dimensional simplex it is det A = 0 and one can write Thus for the relative interior of a simplex it is A −1 1 µ > 0 and every of the n + 1 rows of A −1 1 µ = 0 describes a facet of the simplex.
Vertices of the simplex are the intersection of n facets. Thus, a vertex v of P is the solution of n rows of the linear equation system A −1 x = 0. Clearly the k-th column of A solves the system A −1 x = 0 except of the k-th row. Hence, the intersection of n facets is the k-th column of A, where k is the index which belongs to the facet which is not involved in the intersection. In a simplex that means that the the k-th facet is opposite to the k-th vertex.
In the case of a general convex polytope the connection between vertices and intersection of halfspaces is more complex. Nevertheless, one can establish an inequality which holds for all points in the polytope. Lemma 6.2: Let P = Conv(A) ⊂ R n a full dimensional polytope with N vertices and ν = (ν 0 , ν) ∈ C n+1 , with Re ν 0 > 0 as before. Consider further an index set of vertices σ ⊂ {1, . . . , N } which corresponds to a simplex, det A σ = 0. Then the following equivalence holds where ν/ν 0 := (ν 1 /ν 0 , . . . ν n /ν 0 ) is the componentwise division.
Proof. From the vertex definition of polytopes it is for k,k ∈ R N ≥0 . Since the polytope is full dimensional, there is a triangulation and one can divide the polytope in A = (A σ , Aσ) where A σ has an inverse. Thus the lemma follows by multiplying by the inverse A −1 σ . This lemma establishes a connection between the signs of the matrix A −1 σ Aσ and the geometry of the polytope. Since the signs of A −1 σ Aσ determine the kinematic regions where the series representations (theorem 3.2) converge fast, this lemma can be useful to analyze triangulations in a numerical application.
Furthermore, the matrices A σ and Aσ which appear in the triangulation procedure has some special properties, which also guaranties the convergence of Γ-series. whereσ is the complement of σ and ν = (ν 0 , ν 1 , . . . , ν n ) ∈ C n+1 an arbitrary complex vector.
Proof. The columns of the matrix B = −A −1 σ Aσ 1 ∈ Q N ×r are a basis of the lattice L := ker(A) ∩ Z, which can be verified by direct computation AB = 0. It follows particularly that j B jk = 0 and thereby equation (6.4). The second statement (6.5) follows trivially from j, This lemma also implies that the appearing variables x j = (zσ) j i (z σ ) in the Γ-series are dimensionless (without units). Therefore, the unit of the Γ-series will be only determined by the monomial in front of the power series and is equal to [φ σ,k (ν, z)] = [z i ] −ν 0 .
To finish this short section about polytopes, we introduce another method to triangulate convex polytopes. In section 1.2 there was already a method described how to construct triangulations by generic height vectors. Here, we discuss a combinatorial way to triangulate polytopes. This combinatorial approach also allows the structure of subtriangulations in a natural way.
In order to ensure the convexity of polytopes, we will formulate the concept of visibility in mathematical terms. Consider a face F of a convex polytope P ⊂ R n and an arbitrary point x ∈ relint(F ) in the relative interior of F . The face F is visible from another point p / ∈ F , if the line segment [x, p] intersects P only at x [21].
With the concept of visibility, we can add points to a triangulation and obtains an enlarged triangulation. Let T be a regular triangulation of the convex polytope Conv(A). Then the set T = T ∪ {B ∪ {p} | B ∈ T and B is visible from p} (6.6) is a regular triangulation of the convex polytope Conv(A ∪ p) [21]. In doing so, one can create a so-called placing triangulation by starting with an arbitrary point (the triangulation of one point is still the point) and placing then step by step the other points to the previous triangulation. The order of the added points will determine the triangulation.
In practice there are different algorithms available to calculate triangulations of convex polytopes efficiently, e.g. TOPCOM [59] and polymake [25].

Convergence of Γ-series
The Γ-series were introduced in section 2.1 as formal solutions of the GKZ hypergeometric system. Using standard arguments, we proof, that those Γ-series always have a nonvanishing region of convergence. The proof is roughly orientated towards [71]. In order to show the convergence of the Γ-series one has to estimate the summands. As an application of the Stirling formula one can state the following lemma: Lemma 6.4 [similar to [71]]: For every C ∈ C there are constants κ, R ∈ R >0 independent of M , such that In order to estimate products of "self-exponential" functions a a the following inequality is helpful.
Lemma 6.5: The following estimation holds for a where N ∈ N.
Proof. The left inequality is trivially true. The right inequality will be proven by considering different cases of N . For N = 2 without loss of generality it is a 1 a 2 := r ≥ 1. The latter is then equivalent to the Bernoulli inequality a a 1 1 a a 2 2 ≥ a 1 + a 2 2 a 1 +a 2 ⇔ 1 + r − 1 1 + r r+1 ≥ r . (6.12) For N = 2 n with n ∈ N the lemma can be reduced to the case N = 2 by an iterative use. All the other cases can be reduced to a 2 n -case by adding the mean value ∆ : Now, one can estimate the summands of the Γ-series to find a convergent majorant.
Theorem 6.6 [Convergence of the Γ-series [28]]: There is always a non-vanishing interval I = (a, b) = ∅ with a, b ∈ R N >0 such that the Γ-series of equation ( Due to the form of Λ k it is Cλ ∈ Z N ×r . Furthermore by lemma 6.3 it is also N i=1 C ij = 0, which will imply D := C ij >0 C ij λ j = − C ij >0 C ij λ j = 1 2 | ij C ij λ j |. Thus, one can estimate by lemma 6.4 With lemma 6.5 one can estimate further by splitting the product in the N 1 factors with positive exponent − r j=1 C ij λ j > 0 and the N 2 factors with negative exponent − r j=1 C ij λ j < 0. With R max = max i R i we get our final estimation Thus, the Γ-series can be estimated by a geometric series and there is always a nonvanishing region of absolute convergence.
Therefore, Γ-series have always a non-vanishing region of convergence and provide actual solutions of the GKZ hypergeometric system. In addition to it, it was shown in [26] that all Γ-series of a triangulation have a common region of absolute convergence.

Characteristics of Several Simple Graphs
In order to classify different Feynman graphs with regard to their complexity from a hypergeometric perspective, we collate some characterizing numbers for several standard Feynman graphs in table 1. According to theorem 3.2 one expect vol 0 ∆ G series with dimension r := N − n − 1 in a hypergeometric series representation. Every of those series contains (n + 1) Γ-functions. Further, every regular, unimodular triangulation T u gives a way to write the Feynman integral as Horn hypergeometric functions.
If the number of master integrals C = (−1) n χ((C ) n \ {G = 0}) [6] is less than the volume of the Newton polytope C < vol 0 ∆ G one expects linear dependencies in the limit from the generic Feynman integral to the "ordinary", non-generic Feynman integral according to section 3.2.