Velleda Baldoni
email: baldoni@mat.uniroma2.it

Nicole Berline
email: nicole.berline@math.polytechnique.fr

Jesús A De Loera
email: deloera@math.ucdavis.edu

Matthias Köppe

Michèle Vergne
email: vergne@math.polytechnique.fr

M K Öppe

How to Integrate a Polynomial over a Simplex

de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

Introduction

Let ∆ be a d-dimensional rational simplex inside R n and let f ∈ Q[x 1 , . . . , x n ] be a polynomial with rational coefficients. We consider the problem of how to efficiently compute the exact value of the integral of the polynomial f over ∆, which we denote by ∆ f dm. We use here the integral Lebesgue measure dm on the affine hull ∆ of the simplex ∆, defined below in section 2.1. This normalization of the measure occurs naturally in Euler-Maclaurin formulas for a polytope P , which relate sums over the lattice points of P with certain integrals over the various faces of P . For this measure, the volume of the simplex and every integral of a polynomial function with rational coefficients are rational numbers. Thus the result has a representation in the usual (Turing) model of computation. This is in contrast to other normalizations, such as the induced Euclidean measure, where irrational numbers appear.

The main goals of this article are to discuss the computational complexity of the problem and to provide methods to do the computation that are both theoretically efficient and have reasonable performance in concrete examples.

Computation of integrals of polynomials over polytopes is fundamental for many applications. We already mentioned summation over lattice points of a polytope. They also make an appearance in recent results in optimization problems connected to moment matrices [START_REF] Laurent | Sums of squares, moment matrices, and optimization over polynomials[END_REF]. These integrals are also commonly computed in finite element methods, where the domain is decomposed into cells (typically simplices) via a mesh and complicated functions are approximated by polynomials (see for instance [START_REF] Zienkiewicz | The finite element method[END_REF]). When studying a random univariate polynomial p(x) whose coefficients are independent random variables in certain intervals, the probability distribution for the number of real zeros of p(x) is given as an integral over a polytope [START_REF] Bharucha-Reid | Random polynomials[END_REF]. Integrals over polytopes also play a very important role in statistics, see, for instance, [START_REF] Lin | Marginal likelihood integrals for mixtures of independence models[END_REF]. Remark that among all polytopes, simplices are the fundamental case to consider for integration since any convex polytope can be triangulated into finitely many simplices.

Regarding the computational complexity of our problem, one can ask what happens with integration over arbitrary polytopes. It is very educational to look first at the case when f is the constant polynomial 1, and the answer is simply a volume. It has been proved that already computing the volume of polytopes of varying dimension is #P-hard [START_REF] Dyer | On the complexity of computing the volume of a polyhedron[END_REF][START_REF] Brightwell | Counting linear extensions[END_REF][START_REF] Khachiyan | Complexity of polytope volume computation, New trends in discrete and computational geometry[END_REF][START_REF] Lawrence | Polytope volume computation[END_REF], and that even approximating the volume is hard [START_REF] Elekes | A geometric inequality and the complexity of computing volume[END_REF]. More recently in [START_REF] Rademacher | Approximating the centroid is hard[END_REF] it was proved that computing the centroid of a polytope is #P-hard. In contrast, for a simplex, the volume is given by a determinant, which can be computed in polynomial time. One of the key contributions of this paper is to settle the computational complexity of integrating a non-constant polynomial over a simplex.

Before we can state our results let us understand better the input and output of our computations. Our output will always be the rational number ∆ f dm in the usual binary encoding. The d-dimensional input simplex will be represented by its vertices s 1 , . . . , s d+1 (a Vrepresentation) but note that, in the case of a simplex, one can go from its representation as a system of linear inequalities (an H-representation) to a V -representation in polynomial time, simply by computing the inverse of a matrix.

Thus the encoding size of ∆ is given by the number of vertices, the dimension, and the largest binary encoding size of the coordinates among vertices. Computations with polynomials also require that one specifies concrete data structures for reading the input polynomial and to carry on the calculations. There are several possible choices. One common representation of a polynomial is as a sum of monomial terms with rational coefficients. Some authors assume the representation is dense (polynomials are given by a list of the coefficients of all monomials up to a given total degree r), while other authors assume it is sparse (polynomials are specified by a list of exponent vectors of monomials with non-zero coefficients, together with their coefficients). Another popular representation is by straight-line programs. A straight-line program which encodes a polynomial is, roughly speaking, a program without branches which enables us to evaluate it at any given point (see [START_REF] Bürgisser | Algebraic complexity theory[END_REF][START_REF] Matera | Integration of multivariate rational functions given by straight-line programs[END_REF] and references therein). As we explain in Section 2, general straightline programs are too compact for our purposes, so instead we restrict to a subclass we call single-intermediate-use (division-free) straight-line programs or SIU straight-line programs for short. The precise definition and explanation will appear in Section 2 but for now the reader should think that polynomials are represented as fully parenthesized arithmetic expressions involving binary operators + and ×. Now we are ready to state our first result.

Theorem 1 (Integrating general polynomials over a simplex is hard).

The following problem is NP-hard. Input: (I 1 ) numbers d, n ∈ N in unary encoding, (I 2 ) affinely independent rational vectors s 1 , . . . ,

s d+1 ∈ Q n in binary encoding, (I 3 ) an SIU straight-line program Φ encoding a polynomial f ∈ Q[x 1 , . . . , x n ]
with rational coefficients. Output, in binary encoding:

(O 1 ) the rational number ∆ f dm, where ∆ ⊆ R n is the simplex with vertices s 1 , . . . , s d+1 and dm is the integral Lebesgue measure of the rational affine subspace ∆ .

But we can also prove the following positive results.

Theorem 2 (Efficient integration of polynomials of fixed effective number of variables). For every fixed number D ∈ N, there exists a polynomial-time algorithm for the following problem. Input: (I 1 ) numbers d, n, M ∈ N in unary encoding, (I 2 ) affinely independent rational vectors s 1 , . . . , s d+1 ∈ Q n in binary encoding, (I 3 ) a polynomial f ∈ Q[X 1 , . . . , X D ] represented by either an SIU straight-line program Φ of formal degree at most M, or a sparse or dense monomial representation of total degree at most M, (I 4 ) a rational matrix L with D rows and n columns in binary encoding, the rows of which define D linear forms x → ℓ j , x on R n . Output, in binary encoding:

(O 1 ) the rational number ∆ f ( ℓ 1 , x , . . . , ℓ D , x ) dm, where ∆ ⊆ R n is the simplex with vertices s 1 , . . . , s d+1 and dm is the integral Lebesgue measure of the rational affine subspace ∆ .

In particular, the computation of the integral of a power of one linear form can be done by a polynomial time algorithm. This becomes false already if one considers powers of a quadratic form instead of powers of a linear form. Actually, we prove Theorem 1 by looking at powers Q M of the Motzkin-Straus quadratic form of a graph.

As we will see later, when its degree is fixed, a polynomial has a polynomial size representation in either the SIU straight-line program encoding or the sparse or dense monomial representation and one can switch between the two representations efficiently. The notion of formal degree of an SIU straight-line program will be defined in Section 2.

Corollary 3 (Efficient integration of polynomials of fixed degree). For every fixed number M ∈ N, there exists a polynomial-time algorithm for the following problem. Input:

(I 1 ) numbers d, n ∈ N in unary encoding, (I 2 ) affinely independent rational vectors s 1 , . . . ,

s d+1 ∈ Q n in binary encoding, (I 3 ) a polynomial f ∈ Q[x 1 , . . . ,
x n ] represented by either an SIU straight-line program Φ of formal degree at most M, or a sparse or dense monomial representation of total degree at most M. Output, in binary encoding:

(O 1 ) the rational number ∆ f (x) dm, where ∆ ⊆ R n is the simplex with vertices s 1 , . . . , s d+1 and dm is the integral Lebesgue measure of the rational affine subspace ∆ .

Actually, we give two interesting algorithms that prove Corollary 3. First, we simply observe that a monomial with total degree M involves at most M variables. Our second algorithm is related to the polynomial Waring problem: we decompose a homogeneous polynomial of total degree M into a sum of M-th powers of linear forms.

In [START_REF] Lasserre | The multidimensional version of b a x p dx[END_REF] Lasserre and Avrachenkov compute the integral ∆ f (x) dm when f is a homogeneous polynomial, in terms of the corresponding polarized symmetric multilinear form (Proposition 16). We show that their formula also leads to a proof of Corollary 3. Furthermore, several other methods can be used for integration of polynomials of fixed degree. We discuss them in Section 4.

This paper is organized as follows: After some preparation in Section 2, the main theorems are proved in Section 3. In Section 4, we discuss extensions to other convex polytopes and give a survey of the complexity of other algorithms. Finally, in Section 5, we report on a few computational experiments.

Preliminaries

In this section we prepare for the proofs of the main results.

2.1. Integral Lebesgue measure on a rational affine subspace of R n . On R n itself we consider the standard Lebesgue measure, which gives volume 1 to the fundamental domain of the lattice Z n . Let L be a rational linear subspace of dimension d ≤ n. We normalize the Lebesgue measure on L, so that the volume of the fundamental domain of the intersected lattice L ∩ Z n is 1. Then for any affine subspace L + a parallel to L, we define the integral Lebesgue measure dm by translation. For example, the diagonal of the unit square has length 1 instead of √ 2.

Table 1. The representation of (

x 2 1 + • • • + x 2 n ) k as a straight-line program

Intermediate

Comment

q 1 = 0 q 2 = x 1 q 3 = q 2 • q 2 q 4 = q 1 + q 3
Thus q 4 = x 2 1 . q 5 = x 2 q 6 = q 5 • q 5 q 7 = q 4 + q 6

Now q 7 = x 2 1 + x 2 2 . . . .

q 3n-1 = x n q 3n = q 3n-1 • q 3n-1 q 3n+1 = q 3n-2 + q 3n Now q 3n+1 = x 2 1 + • • • + x 2 n . q 3n+2 = 1 q 3n+3 = q 3n+2 • q 3n+1 q 3n+4 = q 3n+3 • q 3n+1
. . .

q 3n+k+2 = q 3n+k+1 • q 3n+1 Final result.
2.2. Encoding polynomials for integration. We now explain our encoding of polynomials as SIU straight-line programs and justify our use of this encoding. We say that a polynomial f is represented as a (division-free) straight-line program Φ if there is a finite sequence of polynomial functions of Q[x 1 , . . . , x n ], namely q 1 , . . . , q k , the so-called intermediate results, such that each q i is either a variable x 1 , . . . , x n , an element of Q, or either the sum or the product of two preceding polynomials in the sequence and such that q k = f . A straight-line program allows us to describe in polynomial space polynomials which otherwise would need to be described with exponentially many monomial terms. For example, think of of the representation of (x is defined recursively in the obvious way, namely as 0 if q i is a constant of Q, as 1 if q i is a variable x j , as the maximum of the formal degrees of the summands if q i is a sum, and the sum of the formal degrees of the factors if q i is a product. The formal degree of Table 2. The representation of a straight-line program for x 2 k , using iterated squaring Intermediate Comment

q 1 = x q 2 = q 1 • q 1 q 3 = q 2 • q 2
. . .

q k+1 = q k • q k Final result.
the straight-line program Φ is the formal degree of the final result q k . Clearly the formal degree gives an upper bound on the degree of the polynomial represented by it.

A favorite example to illustrate the benefits of a straight-line program encoding is that of the symbolic determinant of an n × n matrix. Its dense representation as monomials has size Θ(n!) but it can be computed in O(n 3 ) operations by Gaussian elimination. See the book [START_REF] Bürgisser | Algebraic complexity theory[END_REF] as a reference for this concept.

From a monomial representation of a polynomial of degree M and n variables it is easy to encode it as a straight-line program: first, by going in increasing degree we can write a straight-line program that generates all monomials of degree at most M in n variables. Then for each of them compute the product of the monomial with its coefficient so the length doubles. Finally successively add each term. This gives a final length bounded above by four times the number of monomials of degree at most M in n variables.

Straight-line programs are quite natural in the context of integration. One would certainly not expand (x 2 1 + • • • + x 2 n ) k to carry on numeric integration when we can easily evaluate it as a function. More importantly, straight-line programs are suitable as an input and output encoding and data structure in certain symbolic algorithms for computations with polynomials, like factoring; see [START_REF] Bürgisser | Algebraic complexity theory[END_REF]. Since straight-line programs can be very compact, the algorithms can handle polynomials whose input and output encodings have an exponential size in a sparse monomial representation.

However, a problem with straight-line programs is that this input encoding can be so compact that the output of many computational questions cannot be written down efficiently in the usual binary encoding. For example, while one can encode the polynomial x 2 k with a straight-line program with only k +1 intermediate results (see Table 2), when we compute the value of x 2 k for x = 2, or the integral

2 0 x 2 k dx = 2 2 k +1 /(2 k + 1)
, the binary encoding of the output has a size of Θ(2 k ). Thus the output, given in binary, turns out to be exponentially bigger Table 3. The representation of a single-intermediateuse straight-line program for x 2 k ; note that the iterated squaring method cannot be used

Intermediate

Comment

q 1 = x q 2 = x q 3 = q 1 • q 2
Now q 3 = x 2 , and q 1 and q 2 cannot be used anymore.

q 4 = x q 5 = q 3 • q 4
Thus q 5 = x 3 . . . .

q 2 k+1 -2 = x q 2 k+1 -1 = q 2 k+1 -3 • q 2 k+1 -2 Final result.
than the input encoding. We remark that the same difficulty arises if we choose a sparse input encoding of the polynomial where not only the coefficients but also the exponent vectors are encoded in binary (rather than the usual unary encoding for the exponent vectors).

This motivates the following variation of the notion of straightline program: We say a (division-free) straight-line program is singleintermediate-use, or SIU for short, if every intermediate result is used only once in the definition of other intermediate results. (However, the variables x 1 , . . . , x n can be used arbitrarily often in the definition of intermediate results.) With this definition, all ways to encode the polynomial x 2 k require at least 2 k multiplications. An example SIU straight-line program is shown in Table 3. Clearly single-intermediateuse straight-line programs are equivalent, in terms of expressiveness and encoding complexity, to fully parenthesized arithmetic expressions using binary operators + and ×.

Efficient computation of truncated product of an arbitrary number of polynomials in a fixed number of variables.

The following result will be used in several situations. Lemma 4. For every fixed number D ∈ N, there exists a polynomial time algorithm for the following problem. Input: a number M in unary encoding, a sequence of k polynomials P j ∈ Q[X 1 , . . . , X D ] of total degree at most M, in dense monomial representation. Output: the product P 1 • • • P k truncated at degree M.

Proof. We start with the product of the first two polynomials. We compute the monomials of degree at most M in this product. This takes O(M 2D ) elementary rational operations, and the maximum encoding length of any coefficient in the product is also polynomial in the input data length. Then we multiply this truncated product with the next polynomial, truncating at degree M, and so on. The total computation takes O(kM 2D ) elementary rational operations.

Proofs of the main results

Our aim is to perform an efficient computation of ∆ f dm where ∆ is a simplex and f a polynomial. We will prove first that this is not possible for f of varying degree under the assumption that P = NP. More precisely, we prove that, under this assumption, an efficient computation of ∆ Q M dm is not possible, where Q is a quadratic form and M is allowed to vary.

In the next subsection we present an algorithm to efficiently compute the integral ∆ f dm in some particular situations, most notably the case of arbitrary powers of linear forms.

3.1.

Hardness for polynomials of non-fixed degree. For the proof of Theorem 1 we need to extend the following well-known result of Motzkin and Straus [START_REF] Motzkin | Maxima for graphs and a new proof of a theorem of Turán[END_REF]. In this section, we denote by ∆ the (n

-1)- dimensional canonical simplex { x ∈ R n : x i ≥ 0, n i=1 x i = 1 }
(i,j)∈E(G) x i x j . Then Q G (x) ∞ = 1 2 (1 -1 ω(G)
). Our first result might be of independent interest as it shows that integrals of polynomials over simplices can carry very interesting combinatorial information. This result builds on the theorem of Motzkin and Straus, using the proof of the well-known relation f ∞ = lim p→∞ f p . Lemma 6. Let G be a graph with n vertices and clique number ω(G). Let Q G (x) be the Motzkin-Straus quadratic form. Then for p ≥ 4(e -1)n 3 ln(32n 2 ), the clique number ω(G) is equal to

1 1-2 Q G p .
To prove Lemma 6 we will first prove the following intermediate result.

Lemma 7. For ε > 0 we have

( Q G ∞ -ε)( ε 4 ) (n-1)/p ≤ Q G p ≤ Q G ∞ .
Proof. The right-hand side inequality follows from the normalization of the measure, as

|Q(x)| ≤ Q G ∞ , for all x ∈ ∆.
In order to obtain the other inequality, we use Hölder's inequality ∆ |f g| dm ≤ f p g q , where q is such that 1 p + 1 q = 1. For any (say) continuous function f on ∆, let us denote by ∆(f, ε) the set { x ∈ ∆ : |f (x)| ≥ f ∞ -ε }, and take for g the characteristic function of ∆(f, ε). We obtain

( f ∞ -ε)(vol ∆(f, ε)) 1/p ≤ f p . (1)
Let a be a point of ∆ where the maximum of Q G is attained. Since

∂Q G ∂x i = (i,j)∈E(G) x j we know that 0 ≤ ∂Q G ∂x i ≤ 1 for x ∈ ∆. Since ∆ is convex, we conclude that for any x ∈ ∆, 0 ≤ Q G (a) -Q G (x) ≤ n i=1 |a i -x i | . Thus ∆(Q G , ε) contains the set C ε = { x ∈ ∆ : n i=1 |a i -x i | < ε}. We claim that vol(C ε ) ≥ ( ε 4 ) n-1
. This claim proves the left inequality of the lemma when we apply it to [START_REF] Alexander | Polynomial interpolation in several variables[END_REF].

Consider the dilated simplex ε/2 1+ε/2 ∆ and the translated set

P ε = a 1+ε/2 + ε/2 1+ε/2 ∆. Clearly P ε is contained in ∆. Moreover, for x ∈ P ε , we have n i=1 |a i -x i | ≤ ε 1+ε/2 ≤ ε, hence P ε is contained in C ε . Since vol(∆) = 1 for the normalized measure, the volume of P ε is equal to ( ε/2 1+ε/2 ) n-1 . Hence vol(P ε ) ≥ (ε/4) n-1
. This finishes the proof. Proof of Lemma 6. In the inequalities of Lemma 7, we substitute the relation

Q G ∞ = 1 2 (1 -1 ω(G)
), given by Motzkin-Straus's theorem (Theorem 5). We obtain

( 1 2 (1 - 1 ω(G) ) -ε)(ε/4) n-1 p ≤ Q G p ≤ 1 2 (1 - 1 ω(G)
).

Let us rewrite these inequalities as

1 1 -2 Q G p ≤ ω(G) ≤ 1 1 - 2 Q G p (ε/4) (n-1)/p -2ε . (2)
We only need to prove that for ε = 1 8n 2 and p ≥ 4(e -1)n 3 ln(32n

2 ) we have 0 ≤ L(p) := 1 1 - 2 Q G p (ε/4) n-1 p -2ε - 1 1 -2 Q G p < 1. (3)
Let us write

δ p = Q G p ( 1 (ε/4) n-1 p -1) = Q G p ((32n 2 ) n-1 p -1).
Thus L(p) in Equation 3 becomes now

L(p) = 1 1 -2 Q G p 1 1 -2 δp+ε 1-2 Q G p -1 . (4) Since Q G p ≤ 1 2 (1 -1 ω(G) ) ≤ 1 2 , we have a bound for δ p 0 ≤ δ p ≤ 1 2 ((32n 2 ) n-1 p -1). Let A = ( 4 ε ) n-1 = (32n 2 ) n-1
. Since we assumed p ≥ 4(e-1)n 3 ln(32n 2 ), we have 0

≤ ln A p < 1, hence 0 ≤ A 1/p -1 < (e -1) ln A p . We obtain 0 ≤ δ p ≤ e -1 2 (n -1) log(32n 2 ) p ≤ 1 8n 2 . Since ω(G) ≤ n, we have 1 -2 Q G p ≥ 1/n. Hence we have 2(δ p + ε) 1 -2 Q G p ≤ 1 2n ≤ 1 2 .
Finally for any number 0 < α < 1/2 we have 1 1-α < 1 + 2α, hence applying this fact to Equation 4 with α = 2 δp+ε

1-2 Q G p we get L(p) < 1 1 -2 Q G p 4(δ + ε) 1 -2 Q G p ≤ 4n 2 (δ p + ε) ≤ 1.
This proves Equation 3 and the lemma.

Proof of Theorem 1. The problem of deciding whether the clique number ω(G) of a graph G is greater than a given number K is a well-known NP-complete problem [START_REF] Garey | Computers and intractibility: A guide to the theory of NP-completeness[END_REF]. From Lemma 6 we see that checking this is the same as checking that for p = 4(e -1)n 3 ln(32n 2 ) the integral part of ∆ (Q G ) p dm is less than K p . Note that the polynomial Q G (x) p is a power of a quadratic form and can be encoded as a SIU straightline program of length O(n

3 log n • |E(G)|).
If the computation of the integral ∆ f dm of a polynomial f could be done in polynomial time in the input size of f , we could then verify the desired inequality in polynomial time as well.

3.2. An extension of a formula of Brion. In this section, we obtain several expressions for the integrals ∆ e ℓ dm and

∆ ℓ M 1 1 • • • ℓ M D D dm
, where ∆ ⊂ R n is a simplex and ℓ, ℓ j are linear forms on R n . The first formula, [START_REF] Berline | Local Euler-Maclaurin formula for polytopes[END_REF] in Lemma 8, is obtained by elementary iterated integration on the standard simplex. It leads to a computation of the integral

∆ ℓ M 1 1 • • • ℓ M D
D dm in terms of the Taylor expansion of a certain analytic function associated to ∆ (Corollary 14), hence to a proof of the complexity result of Theorem 2.

In the case of one linear form ℓ which is regular, we recover in this way the "short formula" of Brion as Corollary 11. This result was first obtained by Brion as a particular case of his theorem on polyhedra [START_REF] Brion | Points entiers dans les polyédres convexes[END_REF].

Lemma 8. Let ∆ be the simplex that is the convex hull of (d+1) affinely independent vertices s 1 , s 2 , . . . , s d+1 in R n , and let ℓ be an arbitrary linear form on R n . Then ∆ e ℓ dm = d! vol(∆, dm)

k∈N d+1 ℓ, s 1 k 1 • • • ℓ, s d+1 k d+1 (|k| + d)! , (5)
where |k| = d+1 j=1 k j . Proof. Using an affine change of variables, it is enough to prove [START_REF] Berline | Local Euler-Maclaurin formula for polytopes[END_REF] when ∆ is the d-dimensional standard simplex ∆ st ⊂ R d defined by

∆ st = x ∈ R d : x i ≥ 0, d i=1 x i ≤ 1 .
The volume of ∆ st is equal to 1 d! . In the case of ∆ st , the vertex s j is the basis vector e j for 1 ≤ j ≤ d and s d+1 = 0. Let ℓ, x = d j=1 a j x j . Then [START_REF] Berline | Local Euler-Maclaurin formula for polytopes[END_REF] becomes

∆st e a 1 x 1 +•••+a d x d dx = k∈N d a k 1 1 • • • a k d d (|k| + d)! .
We prove it by induction on d. For d = 1, we have

1 0 e ax dx = e a -1 a = k≥0 a k (k + 1)! . Let d > 1. We write ∆st e a 1 x 1 +•••+a d x d dx = 1 0 e a d x d x j ≥0 x 1 +•••+x d-1 ≤1-x d e a 1 x 1 +•••+a d-1 x d-1 dx 1 . . . dx d-1 dx d .
By the induction hypothesis and an obvious change of variables, the inner integral is equal to

(1 -x d ) d-1 k∈N d-1 (1 -x d ) |k| a k 1 1 • • • a k d-1 d-1 (|k| + d -1)! .
The result now follows from the relation

1 0 (1 -x) p p! e ax dx = k≥0 a k (k + p + 1)! . Remark 9.
Let us replace ℓ by tℓ in [START_REF] Berline | Local Euler-Maclaurin formula for polytopes[END_REF] and expand in powers of t.

We obtain the following formula.

∆ ℓ M dm = d! vol(∆, dm) M! (M + d)! k∈N d+1 ,|k|=M ℓ, s 1 k 1 • • • ℓ, s d+1 k d+1 . ( 6
)
This relation is a particular case of a result of Lasserre and Avrachenkov, Proposition 16, as we will explain in section 4.3 below.

Theorem 10. Let ∆ be the simplex that is the convex hull of (d + 1) affinely independent vertices s 1 , s 2 , . . . , s d+1 in R n .

M ∈N t M (M + d)! M! ∆ ℓ M dm = d! vol(∆, dm) 1 d+1 j=1 (1 -t ℓ, s j ) . (7)
Proof. We apply Formula (6). Summing up from M = 0 to ∞, we recognize the expansion of the right-hand side of ( 7) into a product of geometric series:

M ∈N t M (M + d)! M! ∆ ℓ M dm = d! vol(∆, dm) M ∈N t M k∈N d+1 |,k|=M ℓ, s 1 k 1 • • • ℓ, s d+1 k d+1 .
Corollary 11 (Brion). Let ∆ be as in the previous theorem. Let ℓ be a linear form which is regular w.r.t. ∆, i.e., ℓ, s i = ℓ, s j for any pair i = j. Then we have the following relations.

∆ ℓ M dm = d! vol(∆, dm) M! (M + d)! d+1 i=1 ℓ, s i M +d j =i ℓ, s i -s j . ( 8
) ∆ e ℓ dm = d! vol(∆, dm) d+1 i=1 e ℓ,s i j =i ℓ, s i -s j . (9)
Proof. We consider the right-hand side of (7) as a rational function of t. The poles t = 1/ ℓ, s i are simple precisely when ℓ is regular. In this case, we obtain (8) by taking the expansion into partial fractions. The second relation follows immediately by expanding e ℓ .

When ℓ is regular, Brion's formula is very short, it is a sum of d + 1 terms. When ℓ is not regular, the expansion of (7) into partial fractions leads to an expression of the integral as a sum of residues. Let K ⊆ {1, . . . , d + 1} be an index set of the different poles t = 1/ ℓ, s k , and for k ∈ K let m k denote the order of the pole, i.e.,

m k = # i ∈ {1, . . . , d + 1} : ℓ, s i = ℓ, s k .
With this notation, we have the following formula.

Corollary 12. ∆ ℓ M dm = d! vol(∆, dm) M! (M + d)! k∈K Res ε=0 (ε + ℓ, s k ) M +d ε m k i∈K i =k (ε + ℓ, s k -s i ) m i (10)
Remark 13. It is worth remarking that Corollaries 11 and 12 can be seen as a particular case of the localization theorem in equivariant cohomology (see for instance [START_REF] Berline | Heat kernels and Dirac operators[END_REF]), although we did not use this fact and instead gave a simple direct calculation. In our situation, the variety is the complex projective space CP d , with action of a d-dimensional torus, such that the image of the moment map is the simplex ∆. Brion's formula corresponds to the generic case of a one-parameter subgroup acting with isolated fixed points. In the degenerate case when the set of fixed points has components of positive dimension, the polar parts in [START_REF] Brion | Points entiers dans les polyédres convexes[END_REF] coincide with the contributions of the components to the localization formula.

Finally, we can extend the formula of Theorem 10 on integrating the power of a linear form to the case of a product of powers of several linear forms. Corollary 14. Let ℓ 1 , . . . , ℓ D be D linear forms on R n . We have the following Taylor expansion:

M∈N D t M 1 1 . . . t M D D (|M| + d)! d! vol(∆, dm) ∆ ℓ M 1 1 • • • ℓ M D D M 1 ! . . . M D ! dm = 1 d+1 i=1 (1 -t 1 ℓ 1 , s i -• • • -t D ℓ D , s i ) . (11)
Proof. Replace tℓ with [START_REF] Bharucha-Reid | Random polynomials[END_REF] and take the expansion

t 1 ℓ 1 + • • • + t D ℓ D in
in powers t M 1 1 • • • t M D D . 3.3.
Polynomial time algorithms for polynomial functions of a fixed number of linear forms. We will present an algorithm which, given a polynomial of the particular form f ( ℓ 1 , x , • • • , ℓ D , x ) where f is a polynomial depending on a fixed number D of variables, and ℓ j , x = L j1 x 1 + • • • + L jn x n , for j = 1, . . . , D, are linear forms on R n , computes its integral on a simplex, in time polynomial on the input data. This algorithm relies on Corollary 14.

Proof of Theorem 2. The number of monomials of degree M in D variables is equal to M +D-1

D-1
. Therefore, when D is fixed, the number of monomials of degree at most M in D variables is O(M D ). When the number of variables D of a straight-line program Φ is fixed, it is possible to compute a sparse or dense representation of the polynomial represented by Φ in polynomial time, by a straight-forward execution of the program. Indeed, all intermediate results can be stored as sparse or dense polynomials with O(M D ) monomials. Since the program Φ is single-intermediate-use, the binary encoding size of all coefficients of the monomials can be bounded polynomially by the input encoding size. Thus it is enough to compute the integral of a monomial,

∆ ℓ 1 , x M 1 • • • ℓ D , x M D dm. ( 12
)
From Theorem 14, it follows that

(|M| + d)! d! vol(∆, dm) ∆ ℓ M 1 1 • • • ℓ M D D M 1 ! • • • M D ! dm is the coefficient of t M 1 1 • • • t M D D
in the Taylor expansion of

1 d+1 i=1 (1 -t 1 ℓ 1 , s i -• • • -t D ℓ D , s i )
.

Since D is fixed, this coefficient can be computed in time polynomial with respect to M and the input data, by multiplying series truncated at degree |M|, as explained in Lemma 4. Finally, vol(∆, dm) needs to be computed. If ∆ = conv{s 1 , . . . , s d+1 } is full-dimensional (d = n), we can do so by computing the determinant of the matrix formed by difference vectors of the vertices:

vol(∆, dm) = 1 n! |det(s 1 -s n+1 , . . . , s n -s n+1 )| .
If ∆ is lower-dimensional, we first compute a basis B ∈ Z n×d of the intersection lattice Λ = lin(∆) ∩ Z n . This can be done in polynomial time by applying an efficient algorithm for computing the Hermite normal form [START_REF] Kannan | Polynomial algorithms for computing the Smith and Hermite normal forms of an integer matrix[END_REF]. Then we express each difference vector v i = s is d+1 ∈ lin(∆) for i = 1, . . . , d using the basis B as

v i = Bv ′ i , where v i ∈ Q d . We obtain vol(∆, dm) = 1 d! |det(v ′ 1 , . . . , v ′ d )|
, thus the volume computation is reduced to the calculation of a determinant. This finishes the proof of Theorem 2.

3.4.

Polynomial time algorithms for polynomials of fixed degree. In the present section, we assume that the total degree of the input polynomial f we wish to integrate is a constant M.

Proof of Corollary 3. First of all, when the formal degree M of a straightline program Φ is fixed, it is possible to compute a sparse or dense representation of the polynomial represented by Φ in polynomial time, by a straight-forward execution of the program. Indeed, all intermediate results can be stored as sparse or dense polynomials with O(n M ) monomials. Since the program Φ is single-intermediate-use, the binary encoding size of all coefficients of the monomials can be bounded polynomially by the input encoding size. Now, the key observation is that a monomial of degree at most M depends effectively on D ≤ M variables x i 1 , . . . ,

x i D , thus it is of the form ℓ M 1 1 • • • ℓ M D D
, where the linear forms ℓ j (x) = x i j are the coordinates that effectively appear in the monomial. Thus, Corollary 3 follows immediately from Theorem 2.

In the following, we give another proof of Corollary 3, based on decompositions of polynomials as sums of powers of linear forms.

Alternative proof of Corollary 3. From Corollaries 11 and 12, we derive another efficient algorithm, as follows. The key idea now is that one can decompose the polynomial f as a sum f := ℓ c ℓ ℓ M j with at most 2 M terms in the sum. To handle the case where f is a monomial, we use the well-known identity

x M 1 1 x M 2 2 • • • x Mn n = 1 |M|! 0≤p i ≤M i (-1) |M|-(p 1 +•••+pn) M 1 p 1 • • • M n p n (p 1 x 1 +• • •+p n x n ) |M| , (13)
where

|M| = M 1 + • • • + M n ≤ M.
The number of terms in the sum is the product (M 1 + 1)

• • • (M n + 1), which is bounded by 2 M 1 +•••+Mn = 2 |M| ≤ 2 M
. Since the number of monomials of degree M in n variables is O(n M ), we have a polynomial time algorithm.

Actually, in the implementation of this method, we group together proportional linear forms, thus we often obtain a smaller number of summands.

The problem of finding a decomposition with the smallest possible number of summands is known as the polynomial Waring problem: What is the smallest integer r(M, n) such that a generic homogeneous polynomial f (x 1 , . . . , x n ) of degree M in n variables is expressible as the sum of r(M, n) M-th powers of linear forms? This problem was solved by Alexander and Hirschowitz (see [START_REF] Alexander | Polynomial interpolation in several variables[END_REF], and [START_REF] Brambilla | On the Alexander-Hirschowitz theorem[END_REF] for an extensive survey), but there is no computational or constructive version of this result. In the extreme case, when the polynomial f happens to be the power of one linear form ℓ, one should certainly avoid applying the above decomposition formula to each of the monomials of f . We remark that, when the degree is fixed, we can decide in polynomial time whether a polynomial f , given in sparse or dense monomial representation, is a power of a linear form ℓ, and, if so, construct such a linear form.

Other algorithms for integration and extensions to other polytopes

We conclude with a discussion of how to extend integration to other polytopes and a review of the complexity of other methods to integrate polynomials over polytopes.

4.1.

A formul of Lasserre-Avrachenkov. Another nice formula is the Lasserre-Avrachenkov formula for the integration of a homogeneous polynomial [START_REF] Lasserre | The multidimensional version of b a x p dx[END_REF] on a simplex. As we explain below, this yields a polynomial-time algorithm for the problem of integrating a polynomial of fixed degree over a polytope in varying dimension, thus providing an alternative proof of Corollary 3.

Proposition 16 ([19]

). Let H be a symmetric multilinear form defined on (R d ) M . Let s 1 , s 1 , . . . , s d+1 be the vertices of a d-dimensional simplex ∆. Then one has ∆ H(x, x, . . . , x)dx = vol(∆)

M +d M 1≤i 1 ≤i 2 ≤•••i M ≤d+1 H(s i 1 , s i 2 , . . . , s i M ). ( 14
)
Remark 17. By reindexing the summation in [START_REF] Elekes | A geometric inequality and the complexity of computing volume[END_REF], as H is symmetric, we obtain

∆ H(x, x, . . . , x)dx = vol(∆) M +d M k 1 +•••+k d+1 =M H(s 1 , . . . , s 1 , . . . , s d+1 , . . . , s d+1 ), (15)
where s 1 is repeated k 1 times, s 2 is repeated k 2 times, etc. When H is of the form H = M i=1 ℓ, x i , for a single linear form ℓ, then (15) coincides with Formula (6) in Remark 9. Now any polynomial f which is homogeneous of degree M can be written as f (x) = H f (x, x, . . . , x) for a unique multilinear form

H f . If f = ℓ M then H f = M i=1 ℓ, x i .
Thus for fixed M the computation of H f can be done by decomposing f into a linear combination of powers of linear forms, as we did in the proof of Corollary 3. Alternatively one can use the well-known polarization formula,

H f (x 1 , . . . , x M ) = 1 2 M M! ε∈{±1} M ε 1 ε 2 • • • ε M f M i=1 ε i x i , (16)
Thus from ( 14) we get the following corollary.

Corollary 18. Let f be a homogeneous polynomial of degree M in d variables, and let s 1 , s 2 , . . . , s d+1 be the vertices of a d-dimensional simplex ∆.

Then ∆ f (y) dy = vol(∆) 2 M M! M +d M 1≤i 1 ≤i 2 ≤•••≤i M ≤d+1 ε∈{±1} M ε 1 ε 2 • • • ε M f M k=1 ε k s i k . ( 17
)
We remark that when we fix the degree M of the homogeneous polynomial f , the length of the polarization formula (thus the length of the second sum in ( 17) is a constant. The length of the first sum in [START_REF] Khachiyan | Complexity of polytope volume computation, New trends in discrete and computational geometry[END_REF] is O(n M ). Thus, for fixed degree in varying dimension, we obtain another polynomial-time algorithm for integrating over a simplex.

4.2.

Traditional conversion of the integral as iterated univariate integrals. Let P ⊆ R d be a full-dimensional polytope and f a polynomial. The traditional method we teach our calculus students to compute multivariate integrals over a bounded region requires them to write the integral P f dm is a sum of sequences of one-dimensional integrals

K j=1 b 1j a 1j b 2j a 2j • • • b dj a dj f dx i 1 dx i 2 . . . dx i d (18)
for which we know the limits of integration a ij , b ij explicitly. The problem of finding the limits of integration and the sum has interesting complexity related to the well-known Fourier-Motzkin elimination method (see Chapter One in [START_REF] Ziegler | Lectures on polytopes[END_REF] for a short introduction). Given a system of linear inequalities Ax ≤ b, describing a polytope P ⊂ R d , Fourier-Motzkin elimination is an algorithm that computes a new larger system of inequalities Âx ≤ b with the property that those inequalities that do not contain the variable x d describe the projection of P into the hyperplane x d = 0. We will not explain the details, but Fourier-Motzkin elimination is quite similar to Gaussian elimination in the sense that the main operations necessary to eliminate the last variable x d require to rearrange, scale, and add rows of the matrix (A, b), but unlike Gaussian elimination, new inequalities are added to the system. It was first observed by Schechter [START_REF] Schechter | Integration over a polyhedron: An application of the Fourier-Motzkin elimination method[END_REF] that Fourier-Motzkin elimination provides a way to generate the traditional iterated integrals. More precisely, let us call P d the projection of P into the the hyperplane x d = 0. Clearly when integrating over a polytopal region we expect that the limits of integration will be affine functions. From the output of Fourier-Motzkin Âx ≤ b, we have x ∈ P if and only if (x 1 , . . . , x d-1 ) ∈ P d and for the first k + r inequalites of the system

x d ≤ bi - d-1 j=1 âij x j = A u i (x 1 , . . . , x d-1 )
for i = 1, . . . k as well as

x d ≥ bk+i - d-1 j=1 âk+ij x j , = A l i (x 1 , . . . , x d-1 )
for i = 1, . . . , r. Then, if we define m(x 1 , . . . , x d ) = max{A l j (x 1 , . . . , x d-1 ), j = 1, . . . , r} and M(x 1 , . . . , x d ) = min{A u j (x 1 , . . . , x d-1 ), j = 1, . . . , r}, we can write

P f (x) dm = P d M m f (x) dx 1 dx 2 • • • dx d
Finally the convex polytope P d can be decomposed into polyhedral regions where the functions m, M become simply affine functions from among the list. Since the integral is additive we get an expression

P f (x) dm = i,j P ij d A u j A l j f (x) dx 1 dx 2 • • • dx d .
Finally by repeating the elimination of variables we recover the full iterated list in [START_REF] Lasserre | Integration on a convex polytope[END_REF]. As it was observed in [START_REF] Schechter | Integration over a polyhedron: An application of the Fourier-Motzkin elimination method[END_REF], this algorithm is unfortunately not efficient because the iterated Fourier-Motzkin elimination procedure can produce exponentially many inequalities for the description of the projection (when the dimension d varies). Thus the number of summands considered can in fact grow exponentially.

4.3.

Two formulas for integral computation. We would like to review two formulas that are nice and could speed up computation in particular cases although they do not seem to yield efficient algorithms just on their own. First, one may reduce the computation of P f dm to integrals on the facets of P , by applying Stokes formula. We must be careful to use a rational Lebesgue measure on each facet. As shown in ( [START_REF] Barvinok | Exponential integrals and sums over convex polyhedra (in Russian)[END_REF]), we have the following result. Theorem 19. Let {F i } i=1,...,m be the set of facets of a full-dimensional polytope P ⊆ R n . For each i, let n i be a rational vector which is transverse to the facet F i and pointing outwards P and let dµ i be the Lebesgue measure on the affine hull of F i which is defined by contracting the standard volume form of R n with n i . Then It is clear that, by considering the expansion of the analytic function P e a,x dx, we can again obtain an analogous result for polynomials. An alternative proof was provided by [START_REF] Lasserre | Integration on a convex polytope[END_REF]. The above theorem, however, does not necessarily reduce the computational burden because, depending on the representation of the polytope, the number of facets can be large and also the facets themselves can be complicated polytopes. Yet, together with our results we obtain the following corollary for two special cases.

Corollary 20. There is a polynomial-time algorithm for the following problem. Input:

(I 1 ) the dimension n ∈ N in unary encoding, (I 2 ) a list of rational vectors in binary encoding, namely (i) either vectors (h 1 , h 1,0 ), . . . , (h m , h m,0 ) ∈ Q n+1 that describe the facet-defining inequalities h i , x ≤ h i,0 of a simplicial full-dimensional rational polytope P , (ii) or vectors s 1 , . . . , s N ∈ Q n that are the vertices of a simple full-dimensional rational polytope P , (I 3 ) a rational vector a ∈ Q n in binary encoding, (I 4 ) an exponent M ∈ N in unary encoding. Output, in binary encoding, (O 1 ) the rational number

P f (x) dm where f (x) = a, x M
and where dm is the standard Lebesgue measure on R n .

Proof. In the case (i) of simplicial polytopes P given by facet-defining inequalities, we can use linear programming to compute in polynomial time a V -representation for each simplex F i that is a facet of P . By applying Theorem 19 with ta in place of a and extracting the coefficient of t M in the Taylor expansion of the analytic function t → I P (ta), we obtain the formula

P a, x M dx = 1 (M + 1) y, a m i=1 y, n i F i a, x M +1 dµ i ,
which holds for all y ∈ R n with y, a = 0. It is known that a suitable y ∈ Q n can be constructed in polynomial time. The integrals on the right-hand side can now be evaluated in polynomial time using Theorem 2.

In the case (ii) of simple polytopes P given by their vertices, we make use of the fact that a variant of Brion's formula [START_REF] Brightwell | Counting linear extensions[END_REF] actually holds for arbitrary rational polytopes. For a simple polytope P , it takes the following form.

P ℓ M dx = M! (M + n)! N i=1 ∆ i ℓ, s i M +n s j ∈N (s i ) ℓ, s i -s j , (19)
where N(s i ) denotes the set of vertices adjacent to s i in P , and ∆ i = |det(s is j ) j∈N (s i ) |. The right-hand side is a sum of rational functions of ℓ, where the denominators cancel out so that the sum is actually polynomial. If ℓ is regular, that is to say ℓ, s is j = 0 for any i and j ∈ N(s i ), then the integral can be computed by [START_REF] Lasserre | The multidimensional version of b a x p dx[END_REF] which is a very short formula. However it becomes difficult to extend the method which we used in the case of a simplex. Instead, we can do a perturbation. In [START_REF] Lasserre | The multidimensional version of b a x p dx[END_REF], we replace ℓ by ℓ + εℓ ′ , where ℓ ′ is such that ℓ + εℓ ′ is regular for ε = 0. The algorithm for choosing ℓ ′ is bounded polynomially. Then we do expansions in powers of ε as explained in Lemma 4.

4.4.

Triangulation of arbitrary polytopes. It is well-known that any convex polytope can be triangulated into finitely many simplices. Thus we can use our result to extend the integration of polynomials over any convex polytope. The complexity of doing it this way will directly depend on the number of simplices in a triangulation. This raises the issue of finding the smallest triangulation possible of a given polytope. Unfortunately this problem was proved to be NP-hard even for fixed dimension three (see [START_REF] De Loera | Triangulations: Structures and algorithms[END_REF]). Thus it is in general not a good idea to spend time finding the smallest triangulation possible. A cautionary remark is that one can naively assume that triangulations help for nonconvex polyhedral regions, while in reality it does not because there exist nonconvex polyhedra that are not triangulable unless one adds new points. Deciding how many new points are necessary is an NPhard problem [START_REF] De Loera | Triangulations: Structures and algorithms[END_REF].

Some computational experiments

We have written Maple programs to perform some initial experiments with the methods described in this paper. The programs are available at [START_REF] Baldoni | Maple programs accompanying the manuscript How to integrate a polynomial over a simplex[END_REF].

Decomposition of polynomials into powers of linear forms.

We have written a Maple program that decomposes a homogeneous polynomial f as a sum f = ℓ c ℓ ℓ M with at most 2 M powers of linear forms ℓ for each monomial, using the formula (13) from the construction in the alternative proof of Corollary 3. 2 Actually, in the implementation of this method, we group together proportional linear forms, thus we often obtain a smaller number of summands. Table 4 shows the number of terms in the resulting formula when f is a dense polynomial, for different values of M (degree) and n (dimension). In Table 5, we compare these numbers with those guaranteed by the Alexander-Hirschowitz theorem. The numbers illustrate that there is a big space for improvement, using a possible constructive solution of the polynomial Waring problem rather than the simple explicit construction of formula (13). 5.2. Integration of a power of a linear form over a simplex. We have written a Maple program that implements the method of Section 3.3 for the efficient integration of a power of one linear form over a simplex, ∆ ℓ M dm. 3 In a computational experiment, for a given dimension n and degree M, we picked random full-dimensional simplices ∆ and random linear forms ℓ and used the Maple program to compute the integral. Table 6 shows the computation times. 4 2 The decomposition is done by the Maple procedure from monome to linear and is illustrated by the program test-decomposition-powerlinform.maple. 3 The integration is done by the Maple procedure without basic simplex integral. 4 All experiments were done with Maple 10 on a Sun Fire V440 with UltraSPARC-IIIi processors running at 1.6 GHz. The computation times are given in CPU seconds. All experiments were subject to a time limit of 600 seconds per example. 5.3. Integration of a monomial over a simplex. Finally, we tested the algorithms discussed in this paper, for the case of the integration of a monomial over a simplex. For our experiments, we considered random simplices ∆ of dimension n. We consider a random monomial x m where m = (m 1 , . . . , m d ) with exponents m i between M and M + 1. For each choice of n and M, we picked 50 combinations of simplices ∆ and exponent vectors m.

First we decompose a given monomial into a sum of powers of linear forms and integrate each summand using the Maple program discussed above. 5 Table 7 shows the minimum, average, and maximum computation times. Second, we have implemented the Lasserre-Avrachenko polarization formula [START_REF] Khachiyan | Complexity of polytope volume computation, New trends in discrete and computational geometry[END_REF] of Corollary 18, as described in Section 4.1 6 . The running times are shown in Table 8.

Conclusions

We discussed various algorithms for the exact integration of polynomials over simplicial regions. Beside their theoretical efficiency, the simple rough experiments we performed clearly demonstrated that these methods are robust enough to attack rather difficult problems. Our investigations opened several doors for further development, which we will present in a forthcoming paper.

First, we have some theoretical issues expanding on our results. As in the case of volumes and the computation of centroids, it is likely that our hardness result, Theorem 1, can be extended into an inapproximability result as those obtained in [START_REF] Rademacher | Approximating the centroid is hard[END_REF]. Another goal is to study other families of polytopes for which exact integration can be done efficiently. Furthermore, we will present a natural extension of the computation of integrals, the efficient computation of the highest degree coefficients of a weighted Ehrhart quasipolynomial of a simplex. Besides the methods of the present article, these last computations are based on the results of [START_REF] Baldoni | Euler-Maclaurin expansion of Barvinok valuations and Ehrhart coefficients of rational polytopes[END_REF] and [START_REF] Berline | Local Euler-Maclaurin formula for polytopes[END_REF]. Second, our intention has been all along to develop algorithms with a good chance of becoming practical and that allow for clear implementation. Thus we have also some practical improvements to discuss. Instead of Corollary 14, we can use Brion's formula [START_REF] Brightwell | Counting linear extensions[END_REF] and a proper use of iterated Laurent expansion to compute ∆ f dm in polynomial time when the degree of f is fixed. We expect that it will speed up the current techniques we have implemented in practice. Finally, in order to develop practical integration software, it appears that our methods should be coupled with fast techniques for decomposing domains into polyhedral regions (e.g. triangulations).

Acknowledgements

A part of the work for this paper was done while the five authors visited Centro di Ricerca Matematica Ennio De Giorgi at the Scuola Normale Superiore of Pisa, Italy. We are grateful for the hospitality of CRM and the support we received. The third author was also supported by NSF grant DMS-0608785. We are grateful for comments from Jean-Bernard Lasserre and Bernd Sturmfels.

  , and we denote by dm the Lebesgue measure on the hyperplane { x ∈ R n : n i=1 x i = 1 }, normalized so that ∆ has volume 1. For a function f on ∆, denote as usual f ∞ = max x∈∆ |f (x)| and f p = ( ∆ |f | p dm) 1/p , for p ≥ 1. Recall that the clique number of a graph G is the largest number of vertices of a complete subgraph of G. Theorem 5 (Motzkin-Straus). Let G be a graph with n vertices and clique number ω(G). Let Q G (x) be the Motzkin-Straus quadratic form 1 2

Theorem 15 .

 15 A generic homogeneous polynomial of degree M in n variables is expressible as the sum of r(M, n) = n+M -1 M n M-th powers of linear forms, with the exception of the cases r(3, 5) = 8, r(4, 3) = 6, r(4, 4) = 10, r(4, 5) = 15, and M = 2, where r(2, n) = n.

F

  i e a,x dµ i for all a ∈ C n and y ∈ R n such that y, a = 0.

Table 2 .

 2 The number of intermediate results of a straight line program is called its length. To keep track of constants we define the size of an intermediate result as one, unless the intermediate result is a constant in which case its size is the binary encoding size of the rational number. The size of a straight-line program is the sum of the sizes of the intermediate results. The formal degree of an intermediate result q i

	2 1 + • • • + x 2 n ) k

Table 4 .

 4 Decomposition of polynomials into powers of linear forms, monomial-by-monomial construction

				Degree M
	n 1	2	5	10	20	30	40	50
	2 2	5	28	145	963	3047 7019 13503
	3 3	12	162 2094 40041 252367
	4 4	22	572 15476
	5 5	35	1546 78786
	8 8	92 13168
	10 10 145 37192
	15 15 330 252813
	20 20 590
	30 30 1335
	40 40 2380
	50 50 3725

Table 5 .

 5 Decomposition of polynomials into powers of linear forms, Alexander-Hirschowitz bounds • 10 6 2.1 • 10 7 2.1 • 10 8 1.3 • 10 9 15 1 15 776 130751 9.3 • 10 7 7.7 • 10 9 2.2 • 10 11 3.2 • 10 12 20 1 20 2126 1.0 • 10 6 3.4 • 10 9 9.4 • 10 11 7.0 • 10 13 2.3 • 10 15 30 1 30 9276 2.1 • 10 7 9.4 • 10 11 2.0 • 10 15 7.9 • 10 17 1.1 • 10 20 40 1 40 27151 2.1 • 10 8 7.0 • 10 13 7.9 • 10 17 1.3 • 10 21 6.7 • 10 23 50 1 50 63251 1.3 • 10 9 2.3 • 10 15 1.1 • 10 20 6.7 • 10 23 1.0 • 10 27

	Degree M

Table 6 .

 6 Integration of powers of linear forms over simplices

				Degree M
	n	2	10	20	50	100	300 1000
	10	0.0	0.0	0.0	0.1	0.0	0.0	0.0
	20	0.1	0.2	0.1	0.2	0.1	0.2	0.2
	50	1.0	1.4	1.4	1.6	1.6	1.6	1.7
	100	5.1	8.4	8.7	9.2	9.5 10.0 11.0
	200 36	71	84	88	97	110	120
	300 150	320	400	470	520	530
	400 500
	1000

Table 7 .

 7 Integration of a monomial by decomposition into a sum of powers of linear forms

	Degree M

Table 8 .

 8 Integration of a monomial using polarization

	Degree M

V. BALDONI, N. BERLINE, J. A. DE LOERA, M. K ÖPPE, AND M. VERGNE

All algorithms are implemented in the file integration.lib and illustrated in files called test-*.maple. The tables in this section have been created using the files tables-*.maple.

This method is implemented in the Maple procedure without from monome to linear integral.

This method is implemented in the Maple procedure integrate with polarization.