A Feynman Integral and its Recurrences and Associators

We determine closed and compact expressions for the epsilon-expansion of certain Gaussian hypergeometric functions expanded around half-integer values by explicitly solving for their recurrence relations. This epsilon-expansion is identified with the normalized solution of the underlying Fuchs system of four regular singular points. We compute its regularized zeta series (giving rise to two independent associators) whose ratio gives the epsilon-expansion at a specific value. Furthermore, we use the well known one-loop massive bubble integral as an example to demonstrate how to obtain all-order epsilon-expansions for Feynman integrals and how to construct representations for Feynman integrals in terms of generalized hypergeometric functions. We use the method of differential equations in combination with the recently established general solution for recurrence relations with non-commutative coefficients.

1 Introduction Scattering amplitudes describe the interactions of physical states and play an important role to determine physical observables measurable at colliders. In perturbation theory at each order in the expansion scattering amplitudes are comprised by a sum over Feynman diagrams with a fixed number of loops. Each individual Feynman diagram is represented by integrals over loop momenta and integrates to functions, which typically depend on the Lorentz-invariant quantities of the external particles like their momenta, masses and scales. These functions are generically neither rational nor algebraic but give rise to a branch cut structure following from unitarity and the fact that virtual particles may go on-shell.
The class of functions describing Feynman integrals are iterated integrals, elliptic functions and perhaps generalizations thereof. To obtain physical results one is interested in their Laurent series expansion (ǫ-expansion) about the integer value of the space-time dimension D (typically D = 4 − 2ǫ). In the parameter space of the underlying higher transcendental functions this gives rise to an expansion w.r.t. small parameter ǫ around some fixed numbers, which may be integer or rational numbers. Expansions around integer values is in general sufficient for the evaluation of loop integrals arising in massless quantum field theories. However, the inclusion of particle masses in loop integrals or the evaluation of phase space integrals may give rise to half-integer values, cf. e.g. [1].
The module of hypergeometric functions [2] is ubiquitous both in computing tree-level string amplitudes [3] and in the evaluation of Feynman diagrams with loops cf. e.g. [4]. Therefore, finding an efficient procedure to determine power series expansion of these functions is an important problem. Their underlying higher order differential equations lead to recurrence relations, which can be solved explicitly by the methods recently proposed in [5]. This procedure gives an explicit solution to the recurrence relation providing for each order in ǫ a closed, compact and analytic expression. This way we get hands on the ǫ-expansions of this large family of functions [5]. For a subclass of the latter the coefficients of their expansions generically represent multiple polylogarithms (MPLs), which give rise to periods of mixed Tate motives in algebraic geometry. Then, expansions around rational numbers p/q naturally yield q-th roots of unity exp(2πip/q) in the arguments of the MPLs [6,7].
Amplitudes in field-theory very often can be described by certain differential equations or systems thereof with a given initial value problem subject to physical conditions [8]. For generic parameters the corresponding differential equations for generalized hypergeometric functions are Fuchsian differential equations with the regular singular points at 0, 1 and ∞. For a specific subclass (of at least 2 F 1 hypergeometric functions to be specified later) at rational values p/q of parameters (real parameters shifted by p/q) after a suitable coordinate transformation the underlying first order differential equations become a Fuchs system with q + 2 regular singular points at 0, exp(2πir/q) and ∞ with r = 1, . . . , q (Schlesinger system). By properly assigning the Lie algebra and monodromy representations of this linear system of differential equations their underlying fundamental solutions can be matched with the ǫ-expansion of specific generalized hypergeometric functions. As a consequence, each order in ǫ is given by some combinations of MPLs and group-like matrix products carrying the information on the parameter. Furthermore, at special values the latter can be given in terms of their underlying regularized zeta series. The latter give rise to q independent associators, which are defined as ratio of two solutions of the specific differential equation. This way one obtains a very elegant way of casting the full ǫ-expansion of certain generalized hypergeometric functions into the form dictated by the underlying Lie algebra structure and the analytic structure of MPLs. Computing higher orders in the ǫ-expansion is then reduced to simple matrix multiplications. In this work we explicitly work out the case q = 2, which is relevant to a Feynman integral to be discussed in this work and comment on the generic case q = 2.
Feynman integrals are classified according to their topologies. A topology includes all integrals, that consist of the same set of propagators but have different powers thereof. There are integration-by-parts (IBP) identities [9], which allow to reduce all integrals of a topology to a set of master integrals (MIs). One way to compute MIs is to apply the method of differential equations [8] (for reviews, see [10]). The idea of this method is to take derivatives of MIs w.r.t. kinematic invariants and masses. The results are combinations of integrals of the same topology, which can again be written in terms of MIs using the IBP reduction. This way one obtains differential equations for the MIs. In practice the goal is to solve these equations in a Laurent expansion around ǫ = 0. If the differential equations take a suitable form the coefficient functions of the ǫ-expansion can be given in a straightforward iterative form [11]. By replacing integrations with integral operators the iterative solutions can be written as recurrence relations with non-commutative coefficients. Solving these recurrences yields the all-order ǫ-expansions for Feynman integrals. This approach allows to give the Laurent expansions as infinite series explicitly in terms of iterated integrals, thereby providing solutions, which are exact in ǫ. This outcome is equivalent to the representation generic to hypergeometric functions and their generalization, with the advantage that the behaviour of the Feynman integral at ǫ = 0 can be extracted directly.
The present work is organized as follows. In section 2 we show how to systematically and efficiently derive ǫ-expansions of hypergeometric function 2 F 1 ǫa,ǫb ǫc+ 1 2 ;z with one half-integer parameter. In section 2.1 we introduce harmonic polylogarithms (HPLs), integral operators and related objects. In section 2.2 we discuss the differential equation satisfied by the hypergeometric function 2 F 1 ǫa,ǫb ǫc+ 1 2 ;z . A given order of its ǫ-expansion can be derived from this differential equation based on computing successively all lower orders. Next, in section 2.3 by transforming the differential equation to recurrence relations we derive the terms in the power series solution for any given order in ǫ. This strategy yields an all-order expansion, i.e. an infinite Laurent series in ǫ explicitly in terms of iterated integrals. In section 2.4 we investigate the underlying Fuchs system involving four regular singular points describing the hypergeometric function 2 F 1 a,b c+ 1 2 ;z with one half-integer parameter. For the latter we derive analytic solutions in terms of hyperlogarithms and group-like matrix elements encoding the information on the parameters a, b and c. One of these solutions will be matched with the relevant ǫ-expansion of the hypergeometric function 2 F 1 ǫa,ǫb ǫc+ 1 2 ;z . As a consequence the coefficients of the ǫ k -order in the power series expansion are given by a set of the MPLs of degree k entering the fundamental solution supplemented by matrix products of k matrices. Furthermore, we construct the two regularized zeta series generic to the underlying Fuchs system. The latter gives rise to two independent associators whose ratio will be related to the hypergeometric function 2 F 1 a,b at the special point z = 1. In section 3 as an example we discuss the ǫ-expansion of Feynman integrals using a massive one-loop integral. The all-order expansion is derived via the method of differential equations in section 3.1 and this result is used to construct a representation in terms of a hypergeometric function in section 3.2. Benefiting from the fact that the all-order result is exact in ǫ, we can give the hypergeometric representation of the Feynman integral not only in D = 4−2ǫ dimensions but also in general dimensions D. Finally, in section 4 we present for the hypergeometric function 2 F 1 a,b 1− p q +c ;z the underlying Fuchs system involving q + 2 regular singular points. We comment on its generic solutions and the underlying q associators.

Hypergeometric Function with half-integer parameters
The Gaussian hypergeometric function 2 F 1 is given by the power series [2] with parameters a, b, c ∈ R and with the Pochhammer (rising factorial) symbol: The series (2.1) converges absolutely at the unit circle |z| = 1 if the parameters meet the following condition: In the sequel we want to investigate hypergeometric functions (2.1) with some of their parameters shifted by 1/2. This gives rise to the following five possibilities and: It has been shown in [12] that the class of four hypergeometric functions (2.3) can algebraically be written in terms of only one, e.g.: If the functions (2.3) (with the same argument z) have in their parameters additional integer shifts, they can be expressed in terms of (2.3) and first derivatives thereof and some polynomials in the parameters a, b, c and z. Therefore, in this section we shall apply our new technique [5] to solve for recurrences to compute the following ǫ-expansion: The coefficient functions u k (z) are expressible in terms of harmonic polylogarithms (HPLs) with rational coefficients [13]. On the other hand, the type (2.4) involves elliptic functions and will not be discussed here. E.g. we have: 2 F 1 1 2 , 1 2 1 ; z = 2 π K( √ z) with K the elliptic function of first kind.

Harmonic polylogarithms, integral operators and the generalized operator product
HPLs of weight w ≥ 2 are defined recursively as [14] H(m 0 , m; y) = In the next subsection we will encounter differential equations satisfied by functions for which we know their boundary conditions at y = 1. Therefore, it is advantageous to introduce the following integral operators with their integrations starting at the point 1. We shall use a shorter notation for products of these operators:

(2.12)
Since HPLs with argument 1 can have divergences, the second terms on the r.h.s. of eqs. (2.11) reveal that some products of integral operators are ill-defined. This happens when J(1) is the rightmost operator of a product acting on a constant. It can be seen in the following sections that these cases never appear in our results. Divergent HPLs also arise for other operator products, however, they can be removed using the underlying product algebra. Consider for example J(1, 1, 0) and its representation in terms HPLs following from eqs. (2.11): According to the identity the divergences in the last two terms of (2.13) cancel each other.
The essential difference between the integral operators (2.9) and the integrations in the definition (2.6) of HPLs is the lower bound. See also eq. (2.37) or Ref. [15] for polylogarithms with lower integration limits not chosen to be zero.
Products of integral operators (2.9) are closely related to the 'symbol' of the underlying function, together with the prescription that the integration path goes from 1 to z. Indeed, a sequence of integral operators immediately determines a 'symbol', which when integrated against this path turns it into the Chen iterated integral representation of the function (see e.g. Ref. [16] for an introduction for physicists).
Furthermore we introduce the generalized operator product as the sum of all the n α=1 j α j 1 , j 2 , . . . , j n (2.16) possible distinct permutations of non-commutative factors c i , each one appearing j i times (i = 1, 2, . . . , n). We will call the non-negative integers j i indices and the factors c i arguments of the generalized operator product (2.15). Some examples are:

(2.17)
A recursive definition and basic properties of the generalized operator product can be found in [5]. The object (2.15) is useful to handle non-commutative quantities, such as the integral operators (2.9).

Differential equations for hypergeometric functions
In this section we describe some of the achievements originally developed in [13] for the calculation of expansions (2.5). Some of the formulas derived here are the foundation of the method we use in the next section to obtain the all-order expansion.
Applying the differential operator θ := z d dz on the series (2.1), it is easy to show, that hypergeometric functions satisfy: Combining these relations yields the differential equation: For the hypergeometric function the second order differential equation (2.19) in z can be written as a system of two first order differential equations in the variable for ϕ 1 (y) and ϕ 2 (y) := y d dy ϕ 1 (y):

(2.22)
Inserting the expansion (2.5) and using, that the resulting differential equations are valid at any order in ǫ, yields iterative differential equations for the coefficient functions u k (y) and v k (y) := y d dy u k (y): (2.23) Eq. (2.1) can be used to determine boundary conditions at z = 0. This point transforms to the point y = 1 under (2.21), so that the boundary conditions are v k (1) = 0 for k ≥ 0 and u k (y = 1) = 0 for k ≥ 1. The system (2.23) can then be solved iteratively for k ≥ 1 via: (2.24) From eqs. (2.1) and (2.2) follow the lowest orders u 0 (z) = 1 and u k (z) = 0 for k < 0. Starting with those we can use eqs. (2.24) to straightforwardly calculate the ǫ-expansion (2.5) order by order, e.g.: We used eqs. (2.11) to give u 2 (y) and u 3 (y) in terms of HPLs. This can be achieved for u 4 (y) and higher order coefficient functions as well.
The iterative computation of Laurent expansions of hypergeometric functions described in this section has first been presented in [13]. It has been applied to similar function, e.g. to generalized hypergeometric functions with integer parameters in [17]. Eqs. (2.24) are the main results of this subsection. They allow to straightforwardly calculate the expansion (2.5) up to any order in ǫ. However, is is not possible to obtain a given order without knowing its lower orders. This problem will be solved in the next section.

Recurrence relations for generalized hypergeometric functions
In this section we present the all-order expansions for the hypergeometric functions (2.5). The idea is to write the differential equations for the coefficient functions as recurrence relations 2 . This is achieved by replacing the derivatives and integrations in the iterative solutions of these differential equations with differential and integral operators, respectively. In [19] such recurrence relations have been used to calculate α ′ -expansions for generalized hypergeometric functions, which enter open string amplitudes. The recurrence relations are homogeneous and linear. Their coefficients, which involve integral-and differential-operators, are non-commutative. The iterative solutions (2.24) of the differential equations (2.23) are actually equivalent to the recurrence relations. With the latter it becomes, however, more obvious how to calculate Laurent expansions order by order.
More importantly, the general solution for this type of recurrence relations has been presented in [5]. Hence the all-order expansion of (2.5) can now systematically be constructed and straightforwardly be given in a closed form. By all-order we mean a representations for the coefficient functions u k , which include k as a variable and therefore hold for all orders in ǫ. In contrast to that, the method of the previous section allows to compute ǫ-expansions order by order starting with u 0 , u 1 , u 2 and so on. In other words the formula for u k , no matter if in the form of a iterative solution to differential equations as discussed in the previous section or as a recurrence relation as presented in the following, is not given in terms of HPLs or similar functions. Instead, coefficient functions of lower orders are included. On the other hand our all-order result obtained from the solution of the recurrence relation gives coefficient functions for all orders explicitly in terms of iterated integrals.
According to the discussion at the beginning of this section a variety of other hypergeometric functions 2 F 1 with parameters, which differ from those in (2.5) by integer or half-integer shifts, can algebraically be expressed as linear combinations of the function (2.5) and derivatives thereof. Thus our result (2.28) allow to construct the all-order expansions for these functions as well. One of these functions appears in section 3 in a representation of a Feynman integral.

Fuchsian system and its associators
The differential equation (2.19)  can be written as a system of two first order differential equations of Fuchsian class in the variable y defined in (2.21) with the singularities y 0 = 0, y 1 = 1, y 2 = −1, the vector and (cf. also (2.22)): (2.32) In (2.31) the 2 × 2 matrices A i are given by: (2.33) The system (2.31) has four regular singular points at y = 0, 1, −1 and y = ∞ and is known as Schlesinger system. Recently, in Ref. [5] we have thoroughly discussed the latter and its generic solution given in terms of hyperlogarithms [20].
Hyperlogarithms are defined recursively from words w built from an alphabet {w 0 , w 1 , . . .} (with w i ≃ A i ) with l + 1 letters: (2.34) Generically, we have

35) with the MPLs:
Li n 1 ,...,n l (z 1 , . . . , z l ) = 0<k l <...<k 1 More generally one considers the iterated integrals [21] x = −1, because they have essential singularities at these points. However, one can construct unique analytic solutions Φ i normalized at x = x i with the asymptotic behaviour x → x i as we arrive at the following expansions 4 and: The auxiliary variable δ x appears when the integration path contains a point at which the integrand is not analytic. Then δ x denotes the branch above and below the real axis: δ x = +1 for x ∈ H + and δ x = −1 for x ∈ H − [22].
One can compare the two solutions Φ 0 and Φ i referring to the two points x = x 0 and x = x i , i = 1, 2, respectively. The quotient of any two such solutions is a constant noncommutative series known as regularized zeta series giving rise to an associator [23]. By analytic continuation the connection matrix between the solution Φ 0 and Φ 1 is independent of x and gives rise to the first (Drinfeld) associator: (2.44) Furthermore, the connection matrix between the solution Φ 0 and Φ 2 is independent of x and gives rise to the second (Drinfeld) associator: The two associators (2.44) and (2.45) are independent and will be constructed below.

Hypergeometric function expansion from associators
In this subsection we relate our generic solutions to (2.31) to ǫ-expansions of hypergeometric functions with half-integer parameters.
For a, b, c + 1 2 − a, c + 1 2 − b and c + 1 2 non-integers a pair of two independent solutions to (2.19) is characterized by their local exponents at the regular singularities 0, 1 and ∞ described by the Riemann scheme  E.g. at z = 0 we can choose the two independent solutions: which for z → 0 behaves as: The set of solutions (2.61) can be related to the normalized solution Φ 1 (y) given in (2.42) subject to the transformation (2.21), i.e.: This is accomplished by introducing the connection matrix C: By construction for y → 1 the solution Φ 1 (y) behaves as: Comparing (2.63) and (2.66) gives the following connection matrix: (2.67) The first column of C does not depend on z nor y, while the non-vanishing second column is not relevant to us. In fact, with (2.65) the form (2.67) of C allows us to relate u 1 of (2.61) to the first matrix element of the fundamental solution (2.42) as which exactly agrees with the orders displayed in (2.25).
Finally, with (2.55) we can express the function (2.30) evaluated at z = 1 in terms of the two independent associators: The r.h.s. of (2.70) can be found in (2.59). Up to words of length four (2.70) yields where we have used and the non-vanishing fourth order contributions: (2.73) It is worth pointing out the compactness and simplicity of the results (2.69) and (2.71). For a given weight k in (2.5) the coefficients u k (z) are given by the MPLs of degree k entering the fundamental solution (2.41) (subject to (2.21)) supplemented by matrix products of k matrices A 0 , A 1 and A 2 given in (2.33). Furthermore, for a given degree k the coefficients u k (1) are given by the terms of degree k of product (2.59) of the two universal associators (2.51) and (2.52). As a consequence, computing higher orders in the expansion (2.5) is reduced to simple matrix multiplications.

A simple Feynman integral
We use the family of one-loop bubble integrals with massive propagators and p 2 = s as a simple example to demonstrate, how to obtain all-order expansions in ǫ (D = 4 − 2ǫ) for Feynman integrals. This is achieved for a MI of the family (3.1) in section 3.1. In section 3.2 we discuss how to obtain the representation for that MI in terms of a hypergeometric function from the all-order expansion.
The integral I(a 1 , a 2 ) is well-known. It is given in [24] for general dimensions in terms of generalized hypergeometric functions. Furthermore, various representations of I(1, 1) are discussed in [25]. This includes all-order expansions, for example in terms of Nielsen polylogarithms, and representations in terms of hypergeometric functions 2 F 1 . Both references also discuss the case with propagators of different masses.
3.1 All-order Laurent expansion in ǫ IBP identities [26] allow to express any integral I(a 1 , a 2 ) in terms of two MIs (cf. Fig. 1). for a = 1. As second MI we consider the finite integral I(2, 1). We use the Ansatz where the function g(x) depends on the ratio of the squared mass m 2 and the kinematic invariant s through the variable: In the region s < 0, which transforms to 0 < x < 1 under eq. (3.4), the function g(x) satisfies the following first order differential equation in x [8]: Instead of solving this equation directly we consider the Laurent expansion around ǫ = 0: Inserting (3.6) into eq. (3.5) and using the fact, that the resulting differential equation is valid at any order in ǫ, yields first order differential equations for the coefficient functions g k (x): (3.7) We use the limit of vanishing external momenta, i.e. s = 0 (which corresponds to x = 1), in I(2, 1) s=0 → I(3, 0) and eq. (3.2) to determine the boundary conditions g k (1) = 0. This allows to straightforwardly solve the differential equations (3.7). We can directly calculate g 0 (x) and give an iterative solution for all other coefficient functions: For example we obtain the next orders in terms of HPLs: Representing the Feynman integral by the lowest orders of the ǫ-expansion is sufficient to describe the behaviour around ǫ = 0. But not all informations on the integral are contained in a finite number of orders. Therefore we are interested in the all-order expression for g k (x): an equation valid for all k, which in contrast to eq. (3.9) is given explicitly in terms of iterated integrals. An expression for g k (x) of this type would allow to give the ǫ-expansion (3.6) as an infinite series thereby providing a closed representation of the complete Feynman integral.
Applying the definitions of the integral operators (2.9) in the iterative solution (3.9) yields the first order recurrence relation: (3.12) With the general solution for this type of relations given in [5] and the initial value (3.8) the all-order expression follows directly. The sum is over non-negative integers j 1 , j 2 which fulfil the condition j 1 + j 2 = k. While eq. (3.9) depends on lower orders of the expansion, eq. (3.13) allows to extract the coefficient function of any order in ǫ directly. For example with k = 2 the condition in the sum over j 1 , j 2 has three solutions, which give the following terms: (3.14) In agreement with eq. (3.11) the sum of all three terms gives g 2 (x).
We can now give a compact expression for g(x) in terms of integral operators: Using the recursive definition of the generalized operator product given in [5], it easy to prove, that (3.15) solves the differential equation (3.5). On the other hand, eq. (3.9) or the equivalent recurrence relation (3.12) only solve the differential eq. (3.7) iteratively.
Alternatively to study a Laurent expansion in ǫ in terms of iterated integrals one can look for a solution exact in ǫ in terms of hypergeometric functions and their generalizations. However, for the latter the behaviour of the integral at ǫ = 0 is not obvious. The all-order representation combines the advantages of both approaches. It is given explicitly in terms of iterated integrals, which allows to straightforwardly determine any order of the ǫ-expansion and therefore the behaviour at ǫ = 0. Moreover, the all-order result contains all informations on the Feynman integral in a compact form like a representation in terms of hypergeometric functions. In the next section we give an application for the all-order expression, which confirms this statement.

Representation in terms of a hypergeometric function
Starting with expression (3.15) it is possible to construct the representation of I(2, 1) in terms of a hypergeometric function, because the all-order expression of the latter is known as well. The type of integral operators which are part of (3.15) also appear in the expansion of the hypergeometric function 2 F 1 , which we have considered in section 2. Combining eqs. (2.5) and (2.28) this expansion can be written as: Let us start by comparing the products of integral operators in the all-order expression (3.15) and (3.16). One difference is the operator J(0) to the left of the generalized operator product in (3.16), which is missing in (3.15). Therefore we shall deal with a hypergeometric function, which is the product of the differential operator y d dy and ( ǫ l 1 +l 2 +l 3 +2m (−1) l 2 +m+1 2 l 1 +l 2 which boil down to (2.3) for p q = 1 2 . Again, additional integer shifts in their parameters can be treated in the same way as for the case p q = 1 2 described in the beginning of section 2. The all order ǫ-expansions of the first three functions involve MPLs of q-th roots of unity: r l = exp 2πil q , l = 1, . . . , q . On the other hand, the expansion of the last function is expressible in terms of MPLs of q-th roots of unity times powers of logarithms ln(1 − z) [7]. Again, according to [7] any of the three functions (4.2) can be expressed in terms of the first one (4.1). Therefore, it is sufficient to study the ǫ-expansion of (4.1).
After introducing the new variable [1,7] y = z z − 1 for the function (4.1) can be written as a system of two first order differential equations [7]: 1 q d dy ϕ 1 (y) = y p−1 1 − y q ϕ 2 (y) , − 1 q d dy ϕ 2 (y) = (a 1 + a 2 ) y q−1 1 − y q + c y ϕ 2 (y) + a 1 a 2 y q−p−1 1 − y q ϕ 1 (y) . (4.11) The system (4.9) has q + 2 singularities at y 0 = 0, y i = r i , i = 1, . . . , q and ∞. Note, that for q = 2 the matrices do not boil down to (2.33) because the transformation (4.4) differs from (2.21). Nevertheless, for this case the system (4.9) represents an equivalent representation, which becomes manifest in its solutions subject to the transformation rules of the HPLs under y → 1−y 1+y . In the generic case q > 2 under consideration with the q + 1 singular points x 0 = 0 and x i = r i , i = 1, . . . , q the alphabet A consists of q + 1 letters A * = {w 0 , w 1 , w 2 , . . . w q } and is directly related to the differential forms dx x , q l=1 dx x − r l appearing in (4.9). Now the corresponding hyperlogarithms (2.34) involve the q-th roots of unity (4.3). Again, a grouplike solution Φ to (4.9) taking values in C A with the alphabet A = {A 0 , A 1 , . . . , A q } can be given as formal weighted sum over iterated integrals (with the weight given by the number of iterated integrations) Φ(x) = w∈A * L w (x) w , (4.12) with the MPLs (2.34). The solution (4.12) can be constructed recursively and built by Picard's iterative methods. Again, one can construct q + 1 unique analytic solutions Φ i normalized at x = x i with the asymptotic behaviour x → x i as Their q regularized zeta series (2.46) corresponding to the ratios Z (x i ) (A 0 , A 1 , . . . , A q ) = Φ i (x) −1 Φ 0 (x) , i = 1, . . . , q (4.14) can be computed along the lines presented above and their lowest orders are given in (2.50). Again, the latter can be related to (4.1) at special values [28].