Daniel Barsky

AND Vicente Mu Ñoz
email: vicente.munoz@uma.es

Ricardo P Érez-Marco
email: ricardo.perez-marco@math.imj-prg.fr

ON THE GENESIS OF BBP FORMULAS

Keywords: 2010 Mathematics Subject Classification. 11K16, 11J99 . pi, log(2), normal numbers, BBP formula

We present a general procedure to generate infinitely many BBP and BBP-like formulas for the simplest transcendental numbers. This provides some insight and a better understanding into their nature. In particular, we can derive the main known BBP formulas for π. We can understand why many of these formulas are rearrangements of each other. We also understand better where some null BBP formulas representing 0 come from. We also explain what is the observed relation between some BBP formulas for log 2 and π, that are obtained by taking real and imaginary parts of a general complex BBP formula. Our methods are elementary, but motivated by transalgebraic considerations, and offer a new way to obtain and to search many new BBP formulas and, conjecturally, to better understand transalgebraic relations between transcendental constants.

Introduction

More than 20 years ago, D.H. Bailey, P. Bowein and S. Plouffe ( [START_REF] Bailey | On the Rapid Computation of Various Polylogarithmic Constants[END_REF]) presented an efficient algorithm to compute deep binary or hexadecimal digits of π without the need to compute the previous ones. Their algorithm is based on a series representation for π given by a formula discovered by S. Plouffe,

π = +∞ k=0 1 16 k 4 8k + 1 - 2 8k + 4 - 1 8k + 5 - 1 8k + 6 . (1)
Formulas of similar form for other transcendental constants were known from long time ago, like the classical formula for log 2, that was known to J. Bernoulli,

log 2 = +∞ k=1 1 2 k 1 k .
The reader can find in [START_REF] Bailey | The computation of previously inaccessible digits of π 2 and Catalan's constant[END_REF] an illustration of the way to extract binary digits from this type of formulae.

Many new formulas of this type, named BBP formulas, have been found for π and other higher transcendental constants in the last decades (see [START_REF] Bailey | A compendium of BBP-type formulas for mathematical constants[END_REF], [START_REF] Tam | BBP-type Formula Database[END_REF]). Plouffe's formula, and others for π, can be derived using integral periods (as in [START_REF] Bailey | On the Rapid Computation of Various Polylogarithmic Constants[END_REF]), or more directly using polylogarithm ladder relations at precise algebraic values (as in [START_REF] Broadhurst | Polylogarithmic ladders, hypergeometric series and the ten millionth digits of ζ(3) and ζ[END_REF]), which can be viewed as generalizations of Machin-Störmer relations (see [START_REF] St Örmer | Sur l'application de la théorie des nombres entiers complexes à la solution en nombres rationnels x 1[END_REF] and [START_REF] St Örmer | Solution complète en nombres entiers de l'équation m arctan 1 x + n arctan 1 y = k π 4[END_REF]) for rational values of the arctangent function, and taking its Taylor series expansions. In particular, we can recover in that way Bellard's formula (see Bellard's webpage [START_REF] Bellard | Pi formulas, algorithms and computations[END_REF]), that seems to be the most efficient one for the purpose of computation of deep binary digits of π (see Remark 1.2),

π = 1 2 6 +∞ k=0 (-1) k 2 10k - 2 5 4k + 1 - 1 4k + 3 + 2 8 10k + 1 - 2 6 10k + 3 - 2 2 10k + 5 - 2 2 10k + 7 + 1 10k + 9 . (2)
Many of these formulas are rearrangements of each other, or related by null BBP formulas that represent 0. The origin of null BBP formulas is somewhat mysterious. Most of the formulas of this sort have been found by extensive computer search over parameter space using the PSLQ algorithm to detect integer relations. So their true origin and nature remainded somewhat mysterious. As the authors of [START_REF] Bailey | On the Rapid Computation of Various Polylogarithmic Constants[END_REF] explain:

We found the identity by a combination of inspired guessing and extensive searching using the PSLQ integer relation algorithm.

and in [START_REF] Bailey | Experimental Mathematics: Examples, methods and implications[END_REF] This formula (1) was found using months of PSLQ computations after corresponding but simpler n-th digit formulas were identified for several other constants, including log [START_REF] Bailey | Experimental Mathematics: Recent Developments and Future Outlook, Mathematics Unlimited -2001 and Beyond[END_REF]. This is likely the first instance in history that a significant new formula for π was discovered by a computer.

We note also the observed mysterious numerical relation of BBP formulas for π and log [START_REF] Bailey | Experimental Mathematics: Recent Developments and Future Outlook, Mathematics Unlimited -2001 and Beyond[END_REF].

For the purpose of computation of all digits of π up to a certain order, there are more efficient formulas given by rapidly convergent series of a modular nature, initiated by S. Ramanujan ([23]), that are at the origin of Chudnoskys' algorithm based on Chudnovskys' formula (see [START_REF] Chudnovsky | Approximation and complex multiplication according to Ramanujan[END_REF])

1 π = 12 +∞ k=0
(-1) k (6k)!(545140134k + 13591409) (3k)!(k!) 3 (640320) 3k+3/2 .

Other methods of algorithmic nature include the Borwein quartic algorithm for π (see [START_REF] Borwein | Pi and the AGM -A Study in Analytic Number Theory and Computational Complexity[END_REF]) that approximately quadruples the number of correct digits with each iteration, and the Borwein nonic algorithm for π that approximately yields nine-times the number of correct digits.

A general BBP formula as defined in [START_REF] Bailey | On the random character of fundamental constant expansions[END_REF] for the constant α is a series of the form

α = Q(d, b, m, A) = +∞ k=0 1 b k m l=1 a l (km + l) d ,
where b, d, m, are integers, b ≥ 2, and A = (a 1 , a 2 , . . . , a m ) is an integer vector. The integer d ≥ 1 is the degree of the formula. The classical BBP formula (1) and Bellard formula (2) are of degree 1. We study in this article formulas of degree 1. It would be interesting to extend the present results to get higher degree formulas. The integer b is called the base of the BBP formula, and digits in base b can be computed efficiently. Particular attention has been given to base b = 2 n formulas, as they are useful in computing binary digits. They are called binary BBP formulas. While there are both base 2 and base 3 BBP formulas for some constants like π 2 (see [START_REF] Borwein | Finding and excluding b-ary Machin-type individual digit formulae[END_REF]), no base 3 formula for π is known.

More generally, we can define BBP-like formulas to be of the general form

α = Q(r 0 , r 1 , d, b, m, A) = r 0 + r 1 +∞ k=0 1 b k m l=1 a l (km + l) d , (3)
where r 0 and r 1 are rational numbers. These more general BBP-like formulas have potentially similar computational applications.

But the interest of these formulas is also theoretical. A normal number in base b ≥ 2 is an irrational number α such that its expansion in base b contains any string of n consecutive digits with frequency b -n . These numbers were introduced in 1909 by É. Borel in an article where he proved that Lebesgue almost every number is normal in any base b ≥ 2 ( [START_REF] Borel | Les probabilités dénombrables et leurs applications arithmétiques[END_REF], and the survey [START_REF] Queff Élec | Old and new results on normality[END_REF]). This result is a direct application of Birkhoff's Ergodic Theorem to the dynamical system given by the transformation T : T → T, multiplication by the base T (x) = bx modulo 1, where T = R/Z. The transformation T preserves the Lebesgue measure which is an ergodic invariant measure. It is not difficult to construct explicit normal numbers, and numbers that are not normal, but there is no known example of "natural" transcendental constant that is normal in every base. It is conjectured that this holds for π and other natural transcendental constants, but this remains an open question. It is not even known if a given digit appears infinitely often in the base 10 expansion for π. We note recent results [START_REF] Bailey | Random generators and normal numbers[END_REF][START_REF] Boros | Irresistible integrals[END_REF] where the normality of certain class of constants has been proved, yet not including π.

An approach to prove normality in base b for any transcendental constant which admits a BBP formula in base b is proposed in [START_REF] Bailey | On the random character of fundamental constant expansions[END_REF]. The criterion, named "Hypothesis A", seems related to Furstenberg's "multiplication by 2 and 3" conjecture (see [START_REF] Lagarias | On the normality of arithmetical constants[END_REF]). Only a very particular class of period-like numbers have BBP formulas (for instante, as mentioned before, π 2 does). It is also natural to investigate the class of numbers with a BBP or BBP-like representation.

The main goal of this article is to present a general procedure to generate the most basic BBP and BBP-like formulas of degree 1 that correspond to the simplest transcendental numbers log p and π. With this new procedure we derive the classical formulas, like Bailey-Borwein-Plouffe or Bellard formulas, and understand better their origin, in particular the origin of null formulas, and the relation of BBP formulas for log 2 and π that correspond to take the real or imaginary parts of the same complex formula. We also understand better the redundancy of rearrangements in these formulas, and the method provides a tool to search for more formulas with a more conceptual approach. Although we do not find new BBP formulas, we recover the most important ones and we believe that the method presented can be further developed to discover new ones. We plan to carry this out in the future.

The procedure to generate BBP formulas is elementary and is motivated by considering the bases for first order asymptotics at infinite of Euler Gamma function and higher Barnes Gamma functions and the transalgebraic considerations that play an important role in [START_REF] Mu Ñoz | Weierstrass constants and exponential periods[END_REF] (see also [START_REF] Mu Ñoz | Unified Treatment of Explicit and Trace Formulas via Poisson-Newton formula[END_REF]). To construct these asymptotic bases, we consider the family iterated integrals of 1 s defined by I 0 (s) = 1 s , and for n ≥ 0,

I n+1 (s) = s 1 I n (u) du = . . . = s 1 un 1 . . . u0 1 1 u 0 du 0 . . . du n-1 du n .
It is elementary to check by induction that

I n (s) = A n (s) log s + B n (s),
where A n , B n ∈ Q[s] are polynomials with rational coefficients, with deg

A n = deg B n = n -1, and
A n (s) = s n-1 (n -1)! ,
we have Theorem 1.1. Let s ∈ C, |s -1| < 1, or |s -1| = 1 and n ≥ 2. We have

I n (s) = +∞ j=0 (1 -s) j+n n! j+n n = s n-1 (n -1)! log s + B n (s), or log s = (1 -s) n ns n-1 +∞ j=0 (1 -s) j j+n j - (n -1)! s n-1 B n (s).
Since B n (s) has rational coefficients, we can take s = 1 -1 b and we get a BBP-like formula for log s. Taking suitable complex values for s, and separating real and imaginary parts, we also obtain BBP and BBP-like formulas for π. We prove that formulas for different values of n provide nonobvious rearrangements of the summations, which in part explains the rich "rearrangement algebra" of BBP formulas.

We recover many formulas with this procedure. For instance, all the formulas of log 2 appearing in Wikipedia [START_REF]Natural logarithm of 2[END_REF] are given in ( 6)- [START_REF] Chudnovsky | Approximation and complex multiplication according to Ramanujan[END_REF]. We also get the following classical formulas:

log 2 = 5 6 - +∞ k=1 1 2 k 1 k - 3 k + 1 + 3 k + 2 - 1 k + 3 , log 2 = 2 3 + +∞ k=1 1 16 k 2 8k + 1 8k + 2 + 1/2 8k + 4 + 1/4 8k + 6 , π = 8 3 + 4 +∞ k=1 1 4k + 1 - 1 4k + 3 , π = ∞ k=0 1 16 k 2 8k + 1 + 2 8k + 2 + 1 8k + 3 - 1/2 8k + 5 - 1/2 8k + 6 - 1/4 8k + 7
.

Also combining our formulas we can get some null formulas representing 0, as for example the following one appearing in [START_REF] Bailey | On the Rapid Computation of Various Polylogarithmic Constants[END_REF] 0 = +∞ k=0

This gives some explanations of the mysteries mentioned before. For example, formulas for log 2 and π are related by taking real and imaginary parts of formulas for complex values for s, for example for s = 1+i 2 . Null formulas can appear when comparing our formulas for different complex values of s taking real or imaginary parts. For example for s = 1/2 and s = 1+i 2 we do get the previous null formula. It is natural to ask if all null BBP formulas of degree 1 can be obtaining combining formulas from Theorem 1.1 for different values of s.

Certainly, we also recover the classical BBP formula (1) and Bellard formula (2).

Remark 1.2. We can measure the efficiency of a BBP-like formula (3) for computing a number α as m/ log b, where m is the number of non-zero coefficients in A = (a 1 , a 2 , . . . , a m ), as this measures the number of non-zero summands for going to the next step in the digit computation. Binary BBP formulas, that is when b = 2, are of special relevance, since they allow to compute α in binary form.

In that case, we can take the logarithm in base 2. The efficiency of ( 1) is 1, whereas the efficiency of (2) is 7/10, a 43% faster.

The techniques of this article extend to other bases of iterated functions that we will discuss in future articles. We hope that our approach can be useful in finding more efficient BBP-formulas for π by more powerful algebraic computer search algorithms.

Acknowledgements. We are very grateful to the anonymous referee that has made a large number of interesting suggestions to improve the exposition. We thank also Tomohiro Yamada for pointing out several corrections. The second author was partially supported by Project MINECO (Spain) PGC2018-095448-B-I00.

Laplace-Hadamard regularization of polar parts

The Laplace-Hadamard regularization is related to work in [START_REF] Mu Ñoz | Unified Treatment of Explicit and Trace Formulas via Poisson-Newton formula[END_REF] and [START_REF] Mu Ñoz | Weierstrass constants and exponential periods[END_REF].

For each n ≥ 0 we define the polynomials P 0 = 0, and for n ≥ 1,

P n (s, t) = n-1 k=0 (1 -s) k k! t k .
We also define the iterated primitives of 1/s defined by I 0 (s) = 1 s , and for n ≥ 0,

I n+1 (s) = s 1 I n (u) du = . . . = s 1 un 1 . . . u0 1 1 u 0 du 0 . . . du n-1 du n .
We call the integrals I n (s) the Laplace-Hadamard regularization or the Laplace-Hadamard transform of 1/t n . The functions I n (s) are holomorphic functions in C-]-∞, 0] and have an isolated singularity at 0 with non-trivial monodromy when n ≥ 1. We have a single integral expression for I n (s) as a Laplace-Hadamard regularization:

Proposition 2.1. For n ≥ 0 and ℜs > 0, or ℜs = 0 and n ≥ 2, we have

I n (s) = (-1) n +∞ 0 1
t n e -st -P n (s, t)e -t dt .

Proof. For n = 0 we have

+∞ 0 e -st dt = 1 s ,
and by induction we get the result integrating on the variable s between 1 and s,

I n+1 (s) = s 1 +∞ 0 1
t n e -ut -P n (u, t)e -t dt du , and using that

s 1 e -ut du = - 1 t (e -st -e t ), s 1 P n (u, t)du = - 1 t (P n+1 (s, t) -1).
Note that we have P n (s, t) → e (1-s)t when n → +∞ uniformly on compact sets, and P n (s, t) is the n-th order jet of e (1-s)t at t = 0. So for t → 0 we have e -st -P n (s, t)e -t = O(t n ).

For n = 0 we get the elementary integral

I 0 (s) = +∞ 0 e -st dt = 1 s .
For n = 1 we get the old Frullani integral ( [START_REF] Boros | Irresistible integrals[END_REF] p.98)

I 1 (s) = - +∞ 0 1 t (e -st -e -t ) dt = log s .
We have

I 1 (s) = log s, I 2 (s) = s log s -(s -1), I 3 (s) = s 2 2 log s - 1 4 (s -1)(3s -1), I 4 (s) = s 3 6 log s - 1 36 (s -1)(11s 2 -7s + 2), I 5 (s) = s 4 24 log s - 1 288 (s -1)(25s 3 -23s 2 + 13s -3).
A simple induction shows Proposition 2.2. We have

I n (s) = A n (s) log s + B n (s),
where

A n , B n ∈ Q[s] are polynomials, with deg A n = deg B n = n -1, and
A n (s) = s n-1 (n -1)! .

✷

Regarding the polynomials B n , the relation I ′ n+1 (s) = I n (s) shows that we have

B ′ n+1 (s) = B n (s) - s n-1 n! . ( 5
)
This equation with the condition B n+1 (1) = 0 determines B n+1 uniquely from B n .

We have a formula for B n (see [START_REF] Medina | Iterated primitives of logarithmic powers[END_REF], where I n (s) = f n-1 (x) with x = s -1, and [START_REF] Mathar | Series of reciprocal powers of k-almost primes[END_REF]):

Proposition 2.3. We have for n ≥ 0,

B n+1 (s) = - 1 n! n k=1 n k (H n -H n-k )(s -1) k ,
where H 0 = 0 and H n = 1 + 1 2 + 1 3 + . . . + 1 n are the Harmonic numbers.

Proof. The formula holds for n = 0 and it satisfies B n+1 (1) = 0 and the recurrence relation:

B ′ n+1 (s) = - 1 n! n k=1 k n k (H n -H n-k )(s -1) k-1 = - 1 (n -1)! n k=1 n -1 k -1 (H n -H n-k )(s -1) k-1 = - 1 (n -1)! n-1 k=1 n -1 k -1 (H n-1 -H n-1-(k-1) )(s -1) k-1 - 1 (n -1)! n k=1 1 n (H n -H n-1 )(s -1) k-1 = B n (s) - 1 n! n k=1 n -1 k -1 (s -1) k-1 = B n (s) - 1 n! ((s -1) + 1) n-1 = B n (s) - s n-1 n! .

Now we prove:

Lemma 2.4. For n ≥ 1,

I n+1 (0) = B n+1 (0) = (-1) n+1 n • n! .
We first establish a useful integral representation for harmonic numbers Lemma 2.5.

H n = +∞ 0 e -t -e -(n+1)t 1 -e -t dt .
Proof. We have

H n = n k=1 1 k = n k=1 +∞ 0 e -kt dt = +∞ 0 e -t -e -(n+1)t 1 -e -t dt .
From this it follows Lemma 2.6.

n k=0 n k (-1) k H n-k = (-1) n+1 n .
Proof.

n k=0 n k (-1) k H n-k = +∞ 0 n k=0 n k (-1) k e -t - n k=0 n k (-1) k e -(n-k)t e -t 1 -e -t dt = +∞ 0 (1 -1) n e -t -(-1 + e -t ) n e -t 1 -e -t dt = (-1) n+1 +∞ 0 (1 -e -t ) n-1 e -t dt = (-1) n+1 1 0 x n-1 dx = (-1) n+1 n .
Now we can prove Lemma 2.4.

Proof of Lemma 2.4. We have

B n+1 (0) = - 1 n! n k=1 n k (H n -H n-k )(-1) k = - 1 n! -H n - n k=1 n k (-1) k H n-k = - 1 n! -H n - (-1) n+1 n -H n = (-1) n+1 n • n! . Corollary 2.7. For n ≥ 2, B ′ n+1 (0) = (-1) n (n -1)(n -1)! .
Proof. From (5) we have

B ′ n+1 (0) = B n (0)
, and the result follows from Lemma 2.4. This is related to the following identity with harmonic numbers: Lemma 2.8. For n ≥ 2, we have

n k=0 k n k (-1) k H n-k = (-1) n-1 n n -1 .
Proof. For n ≥ 2, we have

n k=0 k n k (-1) k a n-k = x d dx (a + x) n x=-1 = -n(a -1) n-1 , therefore n k=0 k n k (-1) k H n-k = +∞ 0 -n (1 -1) n-1 e -t + n (e -t -1) n-1 e -t 1 -e -t dt = (-1) n-1 n +∞ 0 (1 -e -t ) n-2 e -t dt = (-1) n-1 n 1 0
x n-2 dx = (-1) n-1 n n -1 .

Egyptian formulas for rational numbers

We start with the simplest case: an Egyptian formula for rationals. The following is an "infinite Egyptian fraction decomposition" for 1 n .

Proposition 3.1 (Infinite Egyptian fraction decomposition). For n ≥ 2, we have

1 n = +∞ j=1 1 j+n+1 n+1
.

Proof. Notice that from Proposition 2.1 we have

I n (s) = (-1) n +∞ 0 1
t n e -st -P n (s, t)e -t dt, with

P n (s, t) = n-1 k=0 (1 -s) k t k k! , hence P n (0, t) = n-1 k=0 t k k! .
So for n ≥ 2, we can develop and exchange the integral and the summation:

(-1) n I n (0) = +∞ 0 e -t t n e t - n-1 k=0 t k k! dt = +∞ 0 e -t +∞ j=0 t j (j + n)! dt = +∞ j=0 j! (j + n)! = 1 n! +∞ j=0 1 j+n n
. Now we have from Lemma 2.4,

I n (0) = B n (0) = (-1) n (n -1)(n -1)! , thus n n -1 = ∞ j=0 1 j+n n
, and the result follows.

As one referee has pointed out to us, Proposition 3.1 follows also by a telescoping sum over

n j+n+1 n+1 = j + n + 1 j+n+1 n+1 - j + n + 2 j+n+2 n+1
, which is found using Gosper's algorithm. We show here that this formula results from our general approach.

BBP-like formulas for log s

In general we have t n e (1-s)t -P n (s, t) dt

= +∞ 0 e -t t n e (1-s)t - n-1 k=0 (1 -s) k t k k! dt = +∞ 0 e -t t n +∞ j=0 (1 -s) j+n t j+n (j + n)! dt = +∞ 0 e -t +∞ j=0 (1 -s) j+n t j (j + n)! dt = +∞ j=0 (1 -s) j+n j! (j + n)! = 1 n! +∞ j=0 (1 -s) j+n j+n n .
Remark 4.2. The formula in Proposition 4.1 also holds for |s -1| = 1 and s = 0, but the convergence of the sum is only conditional. This can be checked by continuity of both sides making |s -1| → 1.

Now, we have

I n (s) = s n-1 (n -1)! log s + B n (s),
and since B n ∈ Q[s] we get,

I n (s) = (-1) n +∞ j=0 (1 -s) j+n n! j+n n = s n-1 (n -1)! log s + B n (s).
In particular, for s ∈ Q we have

+∞ j=0 (1 -s) j+n n! j+n n ∈ Q ⊕ Q log s .
log s = - (n -1)! s n-1 B n (s) + (-1) n (n -1)! s n-1 +∞ j=0 (1 -s) j+n n! j+n n .
We get a group of formulas for log 2 by specializing at s = 2. We have

log 2 = - (n -1)! 2 n-1 B n (2) + (n -1)! 2 n-1 +∞ j=0
(-1) j (j + 1)(j + 2) . . . (j + n)

.

Using the values B 2 (2) = -1, B 3 (2) = -1 2 , B 4 (2) = -8 9 , B 5 (2) = -131 240 , B 6
(2) = -661 3600 , we get:

log 2 = 1 2 + 1 2 ∞ j=1 (-1) j+1 j(j + 1) (6)
log 2 = 5 8 + 1 2 ∞ j=1 (-1) j+1 j(j + 1)(j + 2) (7) log 2 = 2 3 + 3 4 ∞ j=1 (-1) j+1 j(j + 1)(j + 2)(j + 3) (8) log 2 = 131 192 + 3 2 ∞ j=1
(-1) j+1 j(j + 1)(j + 2)(j + 3)(j + 4) ( 9)

log 2 = 661 960 + 15 4 ∞ j=1 (-1) j+1 j(j + 1)(j + 2)(j + 3)(j + 4)(j + 5) (10)
Specializing at s = 1/2, we get the following formula for n ≥ 1.

log 2 = ∞ j=1 1 2 j j (11) log 2 = 1 - ∞ j=1 1 2 j j(j + 1) (12)
log 2 = 1 2 + 2 ∞ j=1 1 2 j j(j + 1)(j + 2) ( 13
) log 2 = 5 6 -6 ∞ j=1 1 2 j j(j + 1)(j + 2)(j + 3) (14) log 2 = 7 12 + 24 ∞ j=1 1 2 j j(j + 1)(j + 2)(j + 3)(j + 4) (15) log 2 = 47 60 -120 ∞ j=1
(-1) j+1 2 j j(j + 1)(j + 2)(j + 3)(j + 4)(j + 5)

All these formulas appear in [START_REF]Natural logarithm of 2[END_REF].

It is customary to write the formulas above by splitting the denominators into simple fractions. For instance, the fourth formula can be written as

log 2 = 5 6 - +∞ j=1 1 2 j 1 j - 3 j + 1 + 3 j + 2 - 1 j + 3 .
If we group for j = 4k, 4k + 1, 4k + 2, 4k + 3, we get

log 2 = 5 6 - +∞ k=1 1 2 4k 1 4k - 3 4k + 1 + 3 4k + 2 - 1 4k + 3 - +∞ k=0 1 2 4k 1/2 4k + 1 - 3/2 4k + 2 + 3/2 4k + 3 - 1/2 4k + 4 - +∞ k=0 1 2 4k 1/
We can obtain many more binary BBP-like formulas. Specializing at s = 3/2 we get the formula for n ≥ 1,

log(3/2) = - 2 3 n-1 (n -1)!B n (3/2) + 2(n -1)!3 n-1 +∞ j=0
(-1) j 2 j (j + 1)(j + 2) . . . (j + n) . As before, the sum can also be written as

log(3/2) = 65 162 + 1 27 +∞ k=1 (-1) k+1 2 k 1 k - 3 k + 1 + 3 k + 2 - 1 k + 3 .
In general, binary BBP formulas can be obtained from Theorem 4.3 by taking s = 1

± 1 2 N , log 2 N ± 1 2 N = -(n -1)! (1 ± 2 -N ) n-1 B n 1 ± 2 -N + (-1) n (n -1)! (1 ± 2 -N ) n-1 +∞ j=0 1 2 N (j+n) (j + 1)(j + 2) . . . (j + n) .
Formulas of this sort are also obtained by Chamberland [START_REF] Chamberland | Binary BBP-formulae for logarithms and generalized Gaussian-Mersenne primes[END_REF].

The numbers 2 and 2 N ± 1, N ≥ 1, generate a multiplicative subgroup of Q * , and for the elements k in that subgroup, we have binary BBP formulas for log k. The first prime that it is not in this subgroup is k = 23. Note that 2 11 -1 = 23 • 89, but these two primes appear always together in the factor decomposition of 2 N -1 when N is a multiple of 11, and do not appear for other values of N . Also they do not appear at all in 2 N + 1, for any natural number N . This can be checked as follows: first 2 11 ≡ 1 (mod 23), so the order of 2 in Z 23 is 11. In particular it cannot be that 2 N ≡ -1 (mod 23), since otherwise 2 2N ≡ 1 (mod 23), and hence 2N |11, so N |11 and thus 2 N ≡ 1 (mod 23). On the other hand, if 2 N ≡ 1 (mod 23) then N is a multiple of 11, and then 23 •89|(2 11 -1)|(2 N -1).

BBP-like formulas for π

We may use Theorem 4.3 for a complex value of s, then we can get BBP-formulas for log k and also for π separating real and imaginary parts. For n = 1 (using Remark 4.2), we have

log s = (s -1) +∞ j=0 (1 -s) j j+1 1 = - +∞ j=1 (1 -s) j j ,
which is the classical series for log s. Make s = 1 + i. We have 0 < ℜ(1

+ i) = 1 < 2 and log(1 + i) = log √ 2 + i π 4 = 1 2 log 2 + i π 4 , and
log(1 + i) = i +∞ j=0 i j j+1 1
.

Separating real and imaginary part and j = 2k or j = 2k + 1 we get two BBP formulas, one for log 2 and the other one for π:

log 2 = 2 +∞ k=0 (-1) k+1 2k+2 1 = +∞ k=0 (-1) k+1 k + 1 and π = 4 +∞ k=0 (-1) k 2k + 1 .
This last formula is just the first Machin formula for π, related to π 4 = arctan 1 .

For general n ≥ 2, we take s = 1 + i, and we have

log(1 + i) = 1 2 log 2 + i π 4 = - (n -1)! (1 + i) n-1 B n (1 + i) + (-1) n (1 -i) n-1 (n -1)! 2 n-1 ∞ j=0 (-i) j+n (j + 1)(j + 2) . . . (j + n) . Let c n = -ℑ (n -1)! (1 + i) n-1 B n (1 + i) , so that π = 4c n + 4(-1) n (n -1)! 2 n-1 n-1 a=0 n -1 a j≡n+a+1 (2)
(-1) (j+n+a+1)/2 (j + 1)(j + 2) . . . (j + n) .

With this machinery at hand, we recover a number of known formulas.

Proposition 5.1 (Leibniz). We have

π = 8 3 + 4 +∞ k=1 1 4k + 1 - 1 4k + 3 .
Proof. We apply the above to n = 2, where we have that

B 2 (1 + i) = -i and c 2 = -ℑ(-i/(1 + i)) = 1/2, thus π = 2 + 2   +∞ j=0 (-1) j (2j + 2)(2j + 3) + +∞ j=0 (-1) j (2j + 1)(2j + 2)   = 2 + 2 +∞ j=0 (-1) j 1 2j + 2 - 1 2j + 3 + 1 2j + 1 - 1 2j + 2 = 2 + 2 +∞ j=0 (-1) j 1 2j + 1 - 1 2j + 3 = 2 + 2 +∞ k=0 1 4k + 1 - 1 4k + 3 - 1 4k + 3 + 1 4k + 5 = 8 3 + 4 +∞ k=1 1 4k + 1 - 1 4k + 3 .
The original BBP formula from [START_REF] Bailey | On the Rapid Computation of Various Polylogarithmic Constants[END_REF] reads as follows:

Theorem 5.2 (Bailey-Borwein-Plouffe). We have

π = +∞ k=0 1 16 k 4 8k + 1 - 2 8k + 4 - 1 8k + 5 - 1 8k + 6 .
Proof. Take s = 1+i 2 , so log s = -1 2 log 2 + i π 4 . Using the formula for n = 1, we have

log 1 + i 2 = - ∞ j=1 (1 -s) j j = - ∞ j=1 (1 -i) j j2 j .
Taking the imaginary part, and agroupping terms for j = 8k + r, r = 1, 2, . . . , 7, 8, we get

π 4 = - ∞ k=0 1 16 k+1 -8 8k + 1 - 8 8k + 2 - 4 8k + 3 + 2 8k + 5 + 2 8k + 6 + 1 8k + 7 so π = ∞ k=0 1 16 k 2 8k + 1 + 2 8k + 2 + 1 8k + 3 - 1/2 8k + 5 - 1/2 8k + 6 - 1/4 8k + 7 . ( 18
)
Similarly, by taking the real part, we get

- 1 2 log 2 = - 71 210 - ∞ k=1 1 16 k+1 16 8k + 8 8k + 1 - 4 8k + 3 - 4 8k + 4 - 2 8k + 5 + 1 8k + 7 so log 2 = ∞ k=0 1 16 k 2 8k + 1 8k + 1 - 1/2 8k + 3 - 1/2 8k + 4 - 1/4 8k + 5 + 1/8 8k + 7 . (19)
Substracting [START_REF] Mu Ñoz | Unified Treatment of Explicit and Trace Formulas via Poisson-Newton formula[END_REF] and our previous formula [START_REF] Lagarias | On the normality of arithmetical constants[END_REF], we get a null formula

0 = +∞ k=0 1 16 k 1 8k + 1 - 1 8k + 2 - 1/2 8k + 3 - 1 8k + 4 - 1/4 8k + 5 - 1/4 8k + 6 + 1/8 8k + 7 (20)
(note that the term k = 0 gives exactly 1/105 = 71/105 -2/3). Adding [START_REF] Mathar | Series of reciprocal powers of k-almost primes[END_REF] to twice this formula, we get

π = +∞ k=0 1 16 k 4 8k + 1 - 2 8k + 4 - 1 8k + 5 - 1 8k + 6 .
In the proof we have proved and used the following null BBP formula that appears in [START_REF] Bailey | On the Rapid Computation of Various Polylogarithmic Constants[END_REF] :

Proposition 5.3. We have +∞ k=0 1 16 k -8 8k + 1 + 8 8k + 2 + 4 8k + 3 + 8 8k + 4 + 2 8k + 5 + 2 8k + 6 - 1 8k + 7 = 0 (21)
Null BBP formulas are very interesting and useful for rewritting BBP formulas. They are obtained by comparing BBP formulas for the same number at different values of s. Proposition 5.4. We have

∞ k=0 1 2 6k 16 6k + 1 - 24 6k + 2 - 8 6k + 3 - 6 6k + 4 + 1 6k + 5 = 0 .
Proof. We use the formulas

log 3 2 = - ∞ k=1 1 2 k (-1) k k = ∞ k=0 1 2 6k 1/2 6k + 1 + -1/4 6k + 2 + 1/8 6k + 3 - 1/16 6k + 4 + 1/32 6k + 5 - 1/64 6k + 6 , log 3 4 = - ∞ k=1 1 4 k 1 k = - ∞ k=0 1 2 6k 1/4 3k + 1 + 1/16 3k + 2 + 1/64 3k + 3 = - ∞ k=0 1 2 6k 1/2 6k + 2 + 1/8 6k + 4 + 1/32 6k + 6 , log 9 8 = - ∞ k=1 1 8 k (-1) k k = ∞ k=0 1 2 6k 1/8 2k + 1 - 1/64 2k + 2 = ∞ k=0 1 2 6k
3/8 6k + 3 -3/64 6k + 6 .

Adding the first two and substracting the third, we get

∞ k=0 1 2 6k 1/2 6k + 1 - 3/4 6k + 2 - 1/4 6k + 3 - 3/16 6k + 4 + 1/32 6k + 5 = 0 .
and multiplying by 32 we get the result.

Finally, we include a proof of Bellard's formula.

Theorem 5.5 (F. Bellard). We have

π = 1 2 6 +∞ k=0 (-1) k 2 10k - 2 5 4n + 1 - 1 4n + 3 + 2 8 10n + 1 - 2 6 10n + 3 - 2 2 10n + 5 - 2 2 10n + 7 + 1 10n + 9 .
Proof. We use the following factorization

1 + i = 2 + i 2 2 7 + i 8 -1
, and taking imaginary parts π 4 = 2ℑ log(1 + i/2) -ℑ log((7 + i)/8) .

For s = (7 + i)/8 and n = 1, we get

ℑ log((7 + i)/8) = -ℑ ∞ j=1 (1 -i) j j8 j = ∞ k=0 1 2 20k 1/8 8k + 1 + 2/8 2 8k + 2 + 2/8 3 8k + 3 - 4/8 5 8k + 5 - 8/8 6 8k + 6 - 8/8 7 8k + 7 = 1 256 ∞ l=0 (-1) l 2 10l 32 4l + 1 + 8 4l + 2 + 1 4l + 3 , (22)
writing j = 8k + r, r = 1, 2, . . . , 8, and then 2k = l.

Now take s = 1 + i/2 and n = 1, to get

ℑ log(1 + i/2) = -ℑ ∞ j=1 (-i) j j2 j = ∞ k=0 1 2 2k+1 (-1) k 2k + 1 = 1 256 ∞ l=0 (-1) l 2 10l 128 10l + 1 - 32 10l + 3 + 8 10l + 5 - 2 10l + 7 + 1/2 10l + 9 . (23)
We substract twice ( 23) minus [START_REF] Queff Élec | Old and new results on normality[END_REF], and use that 2 8 10l+5 -8 4l+2 = -4 10l+5 . Then we get the result.

On the classical BBP form

As defined in [START_REF] Bailey | On the random character of fundamental constant expansions[END_REF] the classical BBP form is

Q(b, d, m, A) = +∞ k=0 1 b k m l=1 a l (km + l) d ,
where b, d, m are integers and A = (a 1 , a 2 , . . . , a m ) is a vector of integers. The degree is d and the base is b. Let us check that with our formula from Theorem 4.3 we get BBP formulas of degree 1. Lemma 6.1. We have

1 j+n n = n l=1 c l j + l ,
where for l = 1, 2, . . . , n, c l is an integer given by

c l = (-1) l-1 n n -1 l -1 .
Proof. As usual, multiply by j + l and set j = -l to get

c l = n! (n -l)(n -l -1) • • • 2 • 1 • (-1)(-2) • • • (-(l -1)) = (-1) l+1 l n l = (-1) l-1 n n -1 l -1 .
We have a general reorganization Lemma that shows that any sum of BBP form with more than m fractions can be reorganized into one with m terms. Lemma 6.2. We have

+∞ j=0 b -j n i=1 c i (mj + i) d = +∞ k=0 b -k m l=1 a l (km + l) d , with a l = i c i b i-l m .
where the sum extends over indexes l + 1 ≤ i ≤ n such that mj + i = mk + l.

Proof. For k = 0, 1, . . . and l = 1, . . . m group the fractions of the sum modulo m with mj + i = mk + l.

These two Lemma prove that the BBP formulas that we get from Theorem 1.1 are of type Q(1, b, 1, (a 1 )).

We can apply this reorganization to the summation in the formula from Theorem 1.1 and get (regrouping the terms with j = k -l + 1 in the third equality),

c l k + 1 (1 -s) k+1+n-l - 0≤k≤l-2 l≤n c l k + 1 (1 -s) k+1+n-l = - 0≤k≤l-2 l≤n c l (1 -s) k+n+1-l k + 1 + +∞ k=0 a k k + 1 (1 -s) k+1 . ( 24
)
with

a k = n l=1 c l (1 -s) n-l . But we have a k = n l=1 c l (1-s) n-l = n l=1 (-1) l-1 n n -1 l -1 (1-s) n-l = (-1) n-1 n(1-(1-s)) n-1 = (-1) n-1 ns n-1 ,
Hence, we recognize in the last sum of (24) log s, so the formula in Theorem 1.1 for n ≥ 2 is a rearrangement of the formula for n = 1 that is the classical Taylor formula for log s

log s = - ∞ k=0 (1 -s) k+1 k + 1 . ( 25
)
We can use this rearrangement to recover directly the formula for the polynomials B n directly:

m≥1 (1 -s) n+m m(m + 1) • • • (m + n) = (-1) n n! m≥1 (1 -s) m+n n k=0 (-1) k n k 1 m + k = (-1) n n! n k=0 (-1) k n k (1 -s) n-k log(s) - k i=1 (1 -s) i i but A n (s) = (-1) n n! n k=0 (-1) k n k (1 -s) n-k = (1 -(1 -s)) n = (-1) n n! s n and B n (s) = (-1) n n! n k=0 (-1) k n k (1 -s) n-k k i=1 (1 -s) i i
which gives after some rearrangment the expression for B n (s)

Formally, there is no extra content in the formulas for the same parameter s but different integers n ≥ 2. However, these rearrangements are computationally useful, and they are not easy to produce. The iterated integrals I n (s) or Proposition 2.2 gives a systematic method to find a family of such resummations. The expression in terms of combinatorical coefficients in the denominator that arise by the iterated integrals in this type of sums can present sometimes some advantages. Of course one is inmediately reminded (even if it is a formula of higher degree) of the famous Apery sum for ζ(3) starting point of his proof of the irrationality of this number.

Appendix. Location of the zeros of the polynomials B n

The application of the formula in Theorem 4.3 to roots of B n , in particular to real roots, gives BBP-like formulas of a special form. We study the location of the roots of B n and the number of real roots.

To understand the polynomials B n (s), we introduce the polynomials C n (x) of degree n -2, for n ≥ 2, defined by

B n (s) = - 1 (n -1)! (s -1)C n (s -1), (26)
so that by Proposition 2.3

C n (x) = n-2 k=0 n -1 k + 1 (H n-1 -H n-k-2 ) x k . ( 27
)
The zeros of D n are those of C n and an extra zero at x = 0. Now we have two interesting equalities:

D ′ n (x) = n-1 k=1 k n -1 k (H n-1 -H n-k-1 )x k-1 = n-1 k=1 (n -1) n -2 k -1 (H n-1 -H n-k-1 )x k-1 = n-2 k=0 (n -1) n -2 k 1 n -1 + H n-2 -H n-k-2 x k = n-2 k=0 (n -1) n -2 k (H n-2 -H (n-1)-k-1 )x k + n-2 k=0 n -2 k x k =(n -1)D n-1 (x) + (1 + x) n-2 , (28)
and

(1 + x)D ′ n -(n -1)D n = (n -1)(1 + x)D n-1 -(n -1)D n + (1 + x) n-1 = n-2 k=0 (n -1) n -2 k (H n-2 -H n-k-2 )x k + n-2 k=0 (n -1) n -2 k (H n-2 -H n-k-2 )x k+1 - n-1 k=0 (n -1) n -1 k (H n-1 -H n-k-1 )x k + (1 + x) n-1 = n-2 k=0 (n -1) n -2 k H n-2 -H n-k-1 + 1 n -k -1 x k + n-1 k=1 (n -1) n -2 k -1 (H n-2 -H n-k-1 )x k - n-1 k=0 (n -1) n -1 k 1 n -1 + H n-2 -H n-k-1 x k + (1 + x) n-1 = n-2 k=0 n -1 n -k -1 n -2 k x k - n-1 k=0 n -1 k x k + (1 + x) n-1 = n-2 k=0 n -1 k x k - n-1 k=0 n -1 k x k + (1 + x) n-1 =(1 + x) n-1 -x n-1 = Q n (x). ( 29
)
Using these equalities, we can prove the following: Proposition 6.4. For n ≥ 2 even, the polynomial C n (x) has no real roots.

For n ≥ 3 odd, the polynomial C n (x) has only one real root and it lies in the interval ]-1, 0[.

Proof.

We want to prove by induction that:

• For n even, x = 0 is the only (simple) zero of D n (x). And D n (x) < 0 for x < 0 and D n (x) > 0 for x > 0. • For n odd, D n has two zeros, at some x 0 ∈] -1, 0[ and at x = 0. And D n (x) > 0 for

x ∈]-∞, x 0 [ ∪ ]0, ∞[ and D n (x) < 0 for x ∈]x 0 , 0[.

Let n be even. We want to prove that D n (x) has only a zero at x = 0. Note that D n (0) = 0 and D ′ n (0) = 1, so D n is increasing at x = 0. For n even we have Q n (x) > 0 everywhere.

• If x ≤ -1 then D n-1 (x) < 0 by induction hypothesis. By (28) we have D ′ n (x) > 0, so it is increasing there. By Lemma 6.3, D n (-1) < 0 so there are no zeros on ]-∞, -1].

• If x > 0 then D n-1 (x) > 0 by induction hypothesis. By (28) we have D ′ n (x) > 0, so it is increasing there. As D n (0) = 0, there are no zeros on ]0, ∞[.

• If x ∈]-1, 0] then Q n (x) > 0. If D n (x) = 0 then (29) says that (1 + x)D ′ n (x) > 0.
So D n is increasing at every zero. As x = 0 is a zero, then this implies that there is only one zero of D n . Now let n be odd. We want to prove that D n (x) has a zero at some x 0 ∈]-1, 0[ and at x = 0, it is positive on ]-∞, x 0 [ ∪ ]0, ∞[ and negative at ]x 0 , 0[. Note that D n (0) = 0 and D ′ n (0) = 1, so it is increasing at x = 0. Note that for n odd we have Q n (x) > 0 for x > -1/2, and Q n (x) < 0 for x < -1/2.

• If x ≤ -1 then D n-1 (x) < 0 by induction hypothesis. By (28) we have D ′ n (x) < 0, so it is decreasing there. By Lemma 6.3 D n (-1) > 0, so there are no zeros on ]-∞, -1].

• If x > 0 then D n-1 (x) > 0 by induction hypothesis. By (28) we have D ′ n (x) > 0, so it is increasing there. As D n (0) = 0, there are no zeros on ]0, ∞[. • If x ∈]-1/2, 0[ then (29) says that (1 + x)D ′ n (x) > (n -1)D n (x). So if there is a zero, D n is increasing. As the last zero before x = 0 cannot be increasing, this last zero has to be x 0 ≤ -1/2. • For x = -1/2, if it was a zero of D n , then it is also a zero of D ′ n because of (29). Then we write x = -1/2 + h, and develop (29) to see that D ′ n (-1/2 + h) > 0 for h > 0 small. But this implies that there must be another zero of D n in ]-1/2, 0[ with decreasing slope, which contradicts the previous item.

• If x ∈]-1, -1/2[ then (29) says that (1 + x)D ′ n (x) < (n -1)D n (x). So if there is a zero, D n is decreasing. There must be at least one zero, but there cannot be two zeros, since there cannot be two decreasing consecutive zeros.

It is relevant to locate the complex zeros of C n (x). The polynomial C 4 has a pair of conjugate complex roots x ≈ -0.68182 ± 0.28386i. The polynomial C 5 has one real root x 0 ≈ -0.61852 and a pair of conjugate complex roots x ≈ -0.73074 ± 0.49200i. The polynomial C 6 has 2 pairs of conjugate complex roots: x ≈ -0.18252 ± 0.39103i, x ≈ -1.2226 ± 0.9258i. We may expect that all roots of C n (x) have ℜx ∈]-∞, 0[. To locate the complex roots of C n (x), we rewrite the differential equation (29) as

(1 + x) -(n-1) D n (x) ′ = (1 + x) -n Q n (x).
Take f n (x) = (1 + x) -(n-1) D n (x), hence df n = (1 + x) -n (1 + x) n-1 -x n-1 dx. We make the change of variables w = x 1+x to get df n = 1-w n-1 1-w dw = (1 + w + .... + w n-2 )dw, and integrating

f n = w + 1 2 w 2 + .... + 1 n -1 w n-1 ,
where we have used that for w = 0, it is x = 0 and hence f n = 0. Note that f n (w) is the truncation of the series -log(1 -w), which is convergent on |w| < 1.

Proposition 6.5. The polynomial f n (w) has no roots in |w| ≤ 1 except w = 0.

Proof. We will look at the polynomial

Q(w) = f n (w)(1 -w)/w = 1 - n-2 k=1 1 k(k + 1) w k - 1 n -1 w n-1 ,
for which we want to check that the only root in the disc |w| ≤ 1 is w = 1. For |w| ≤ 1, we have

n-2 k=2 1 k(k + 1) w k + 1 n -1 w n-1 ≤ n-2 k=2 1 k(k + 1) + 1 n -1 w n-1 = 1 2 .
Then if Q(w) = 0, we have Undoing the change of variables w = x 1+x , we get that all roots of C n (x) are in ℜx < -1 2 . Therefore, with [START_REF] Tam | BBP-type Formula Database[END_REF] we get that the roots of B n (s) are s = 1 and the others lie in ℜs < 1 2 .

Proposition 4 . 1 . 1 t

 411 For |s -1| < 1, or |s -1| = 1 and n ≥ 2, we have I n (s) = (The condition ℜs > 0 ensures the convergence of the integrals and |s -1| < 1, or |s -1| = 1 and n ≥ 2 ensures the convergence of the series, (-1) n I n (s) = +∞ 0 n e -st -P n (s, t)e -t dt =

Theorem 4 . 3 .

 43 Let |s -1| < 1, or |s -1| = 1 and n ≥ 2. Then we have

log 2 = 1 2

 21 -(n -1)! 2 n-1 B n (1/2) + (-1) n (n -1)! +∞ j=0 j+1 (j + 1)(j + 2) . . . (j + n) . Using the values B 2 (1/2) = -2, B 3 (1/2) = -1, B 4 (1/2) = -40 36 , B 5 (1/2) = 7 18 , B 6 (1/2) = -47 225 , we get the formulas:

For instance, n = 4

 4

  |w -2| ≤ 1. Combined with |w| ≤ 1, we have w = 1.

We list the polynomials:

We want to locate the zeros of C n (x).

Lemma 6.3. We have

Proof. The value at x = -1 follows from Lemma 2.4. The value at x = 0 by [START_REF]Natural logarithm of 2[END_REF].

Let


Didier Henrion

Simone Naldi

Mohab Safey

El Din

Algebraic certificates for the truncated moment problem

Keywords: moments, sums of squares, semialgebraic sets, real algebraic geometry, algorithms, complexity

The truncated moment problem consists of determining whether a given finitedimensional vector of real numbers y is obtained by integrating a basis of the vector space of polynomials of bounded degree with respect to a non-negative measure on a given set K of a finite-dimensional Euclidean space. This problem has plenty of applications e.g. in optimization, control theory and statistics. When K is a compact semialgebraic set, the duality between the cone of moments of non-negative measures on K and the cone of non-negative polynomials on K yields an alternative: either y is a moment vector, or y is not a moment vector, in which case there exists a polynomial strictly positive on K making a linear functional depending on y vanish. Such a polynomial is an algebraic certificate of moment unrepresentability. We study the complexity of computing such a certificate using computer algebra algorithms.

Introduction

Problem statement Let x = (x 1 , . . . , x n ) be variables, R[x] be the ring of n-variate real polynomials and for d ∈ N, let R[x] ≤d be the vector space of real polynomials of degree at most d. The multivariate monomial with exponent α = (α 1 , . . . , α n ) ∈ N n is denoted by x α := x α 1 1 • • • x αn n , and its total degree by |α| = α 1 + • • • + α n . For g = (g 1 , . . . , g k ) ∈ R[x] k , the basic semialgebraic set associated with g is S(g ) = {a ∈ R n : g 1 (a) ≥ 0, . . . , g k (a) ≥ 0}.

(1)

Given n, d ∈ N and a sequence of real numbers y = (y α ) α∈N n d indexed by N n d = {α ∈ N n : n i=1 α i ≤ d}, the truncated moment problem (below TMP) is the question of deciding 1 CNRS, LAAS, Université de Toulouse, France. 2 Faculty of Electrical Engineering, Czech Technical University in Prague, Czechia. 3 Université de Limoges, CNRS, XLIM, Limoges, France. 4 Sorbonne Université, LIP6, Paris, France.

whether there exists a nonnegative Borel measure µ on R n , with support in K = S(g ), and such that

y α = K x α d µ, for all α ∈ N n d , (2)
If this is the case, one says that y is moment-representable on K. More generally, the monomial basis can be replaced by another linear basis of R[x] ≤d , (e.g. Chebyshev polynomials).

The TMP is the truncated version of the classical full moment problem [START_REF] Kuhlmann | Positivity, sums of squares and the multi-dimensional moment problem[END_REF][START_REF] Schmüdgen | The moment problem[END_REF].

Overview and state of the art

The TMP is of central importance in data science. It is at the heart of several questions in optimization, control theory or statistics, to mention just a few application domains. It is a key ingredient to the moment-SOS (sum of squares) approach [START_REF] Henrion | The moment-SOS hierarchy[END_REF] which consists of solving numerically non-convex non-linear problems at the price of solving a family of finite-dimensional conic optimization (typically semidefinite programming) problems. Mathematical foundations of the moment problem were recently surveyed in [START_REF] Schmüdgen | The moment problem[END_REF][START_REF] Fialkow | The truncated K-moment problem: a survey[END_REF].

The TMP can be interpreted as a decision problem of the first order theory of the reals, in which case, the input-output data structure is as follows. The input is encoded by a finitedimensional real vector y, whose coordinates are indexed by a basis of the vector space of n-variate polynomials of degree ≤ d, together with finitely-many polynomial inequalities defining a basic semialgebraic set K ⊂ R n . The output is a decision yes/no whether y is obtained by integrating the basis with respect to a non-negative measure supported on K.

When K is compact, the TMP is dual, in the convex analysis sense, to the problem of determining whether a polynomial is non-negative on the semialgebraic set K. This latter problem is at the heart of the development of real algebra in the twentieth century. Whereas deciding nonnegativity of a polynomial of degree at least four can be challenging (NPcomplete problems can be cast as instances of such positivity problems [START_REF] Laurent | Sums of squares, moment matrices and optimization over polynomials[END_REF]), there exist algebraic certificates based on SOS for deciding strict positivity under compactness-like assumptions that can be computed by solving a semidefinite programming (SDP) problem [START_REF] Powers | An algorithm for sums of squares of real polynomials[END_REF]. For positivity over compact basic semialgebraic sets, these certificates have the form of linear combinations with SOS coefficients of polynomials that are explicitly nonnegative on K [START_REF] Putinar | Positive polynomials on compact semi-algebraic sets[END_REF][START_REF] Schmüdgen | The K-moment problem on compact semi-algebraic sets[END_REF]. The size of the SDP is determined by the degree of the SOS-representation.

Seminal papers [START_REF] Ramana | An exact duality theory for semidefinite programming and its complexity implications[END_REF][START_REF] Porkolab | On the complexity of semidefinite programs[END_REF] have started to investigate the complexity of SDP in the context of rational arithmetic. More recent work is based on the determinantal structure of semidefinite programs [START_REF] Henrion | Exact algorithms for linear matrix inequalities[END_REF].

The specific case study of computation of SOS certificates have recently received a lot of attention from the computer algebra community [START_REF] Kaltofen | Exact certification of global optimality of approximate factorizations via rationalizing sums-of-squares with floating point scalars[END_REF][START_REF] Kaltofen | Exact certification in global polynomial optimization via sums-of-squares of rational functions with rational coefficients[END_REF][START_REF] Henrion | Exact algorithms for linear matrix inequalities[END_REF][START_REF] Henrion | SPECTRA-a Maple library for solving linear matrix inequalities in exact arithmetic[END_REF][START_REF] Magron | On exact Reznick, Hilbert-Artin and Putinar's representations[END_REF] especially for the question whether certificates exist over the rational numbers [START_REF] Peyrl | Computing sum of squares decompositions with rational coefficients[END_REF][START_REF] El Din | Computing rational points in convex semialgebraic sets and sum of squares decompositions[END_REF][START_REF] Guo | Computing rational solutions of linear matrix inequalities[END_REF][START_REF] Scheiderer | Sums of squares of polynomials with rational coefficients[END_REF]. In this work we make use of quantifier elimination for determining bounds on the complexity of computing certificates for unrepresentability by using quantitative results from [START_REF] Basu | On the combinatorial and algebraic complexity of quantifier elimination[END_REF]. These also rely on recent advances in the complexity analysis of Putinar Positivstellensatz [START_REF] Baldi | On the effective Putinar's Positivstellensatz and moment approximation[END_REF].

Note that these non-negativity problems can be solved with computer algebra algorithms which are root-finding algorithms, hence which do not provide algebraic certificates of nonnegativity but do provide witnesses (real points) of negativity whenever they exist. The first family of algorithms for doing so is based on the so-called Cylindrical Algebraic Decompo-sition [START_REF] Collins | Quantifier elimination for real closed fields by cylindrical algebraic decomposition[END_REF]. It has complexity which is doubly exponential in the number of variables, which would lead, in the context of TMP, to complexity bounds that are doubly exponential in n + n+d d . The second family of algorithms, named critical point method, initiated by [START_REF] Grigoriev | Solving systems of polynomials inequalities in subexponential time[END_REF], has complexity which is singly exponential in the number of variables (see [START_REF] Basu | Algorithms in real algebraic geometry[END_REF] and references therein).

The purpose of this communication is to leverage on these achievements and initiate the study and development of computer algebra algorithms for solving the truncated moment problem for basic semialgebraic sets.

Overview of the contribution In the wake of the mentioned duality with moments, the existence of SOS-certificates for nonnegative polynomials in the dual side, suggests that similar certificates might be used for the TMP on the primal side.

On the one hand, when a measure exists whose partial moments coincide with y, the measure itself is the natural algebraic proof that allows the user to verify directly that y is moment-representable. On the other hand, this paper shows the existence of explicit algebraic certificates of unrepresentability: these have the form of positive polynomials on K admitting a positivity certificate and orthogonal to the vector y.

Our contribution is based on the fact that the TMP, as a decision problem, is equivalent to the feasibility of a convex conic program in a finite dimensional vector space. More precisely, the question is whether the interior of P(K) d , the cone of polynomials nonnegative on K of degree at most d, intersects the vanishing locus of the Riesz functional L y : R[x] ≤d → R defined by L y ( p α x α ) = p α y α .

When K is compact, Tchakaloff's Theorem [START_REF] Tchakaloff | Formules de cubatures mécaniques à coefficients non négatifs[END_REF] states that y is moment-representable whenever L y is nonnegative on P(K) d , in other words, if the mentioned conic program is weakly feasible. In this case there exists an atomic measure µ = s i=1 c i δ x i whose moment sequence of degree ≤ d is y: such measure is a (real) solution of a highly structured polynomial system of type multivariate Vandermonde, which we do not investigate here. On the other side of the coin, y is not moment-representable exactly when the conic program is strongly feasible: in algebraic terms, this means that there exists a polynomial p ∈ P(K) d , (strictly) positive on K, in the kernel of L y .

In our contribution we study algorithmic aspects of the computation of such unrepresentability algebraic certificates when y is not moment-representable. First, we show that if the quadratic module corresponding to the description of K is archimedean, such certificates exist. We define an integer invariant called the unrepresentability degree which measures the complexity of computing such certificate. We give bounds on such degree that only depend on the input size of our algorithm. When the input vector y is defined over Q, and if it is not moment-representable, we show that there exists a rational certificate of unrepresentability.

Preliminaries

Nonnegative polynomials

Let K ⊂ R n . A polynomial f ∈ R[x] is called nonnegative on K if f (a)
≥ 0 for all a ∈ K and positive if f (a) > 0 for all a ∈ K. We denote by

P(K) d = {f ∈ R[x] ≤d : f (a) ≥ 0, ∀a ∈ K}
the convex cone of polynomials of degree ≤ d, nonnegative on K. If K is semialgebraic, then P(K) d is also semialgebraic by the theorem of Tarski on quantifier elimination over the reals [START_REF] Tarski | A decision method for elementary algebra and geometry[END_REF].

We denote by Σ n ⊂ R[x] the cone of sums of squares of polynomials, and by

Σ n,2d = Σ n ∩ R[x] ≤2d its degree-2d part (remark that Σ n,2d+1 = Σ n,2d ). The cone Σ n,2d is full- dimensional in R[x] ≤2d and is contained in P(R n ) 2d .
Testing membership in cones of nonnegative polynomials over semialgebraic sets can be challenging. Indeed, testing nonnegativity of polynomials of degree ≥ 4 is NP-hard [START_REF] Laurent | Sums of squares, moment matrices and optimization over polynomials[END_REF]Sec. 1.1]. Nevertheless P(K) d contains subcones that can be represented via linear matrix inequalities, thus testing membership in these subcones can be cast as a semidefinite programming (SDP) problem.

Examples of such subsets are quadratic modules: for K = S(g ) as in (1) and d ∈ N, and denoting g 0 := 1, we define the quadratic module associated with g and its truncation of order d respectively by:

Q(g ) = k i=0 σ i g i : σ i ∈ Σ n and Q(g )[d] = k i=0 σ i g i : σ i ∈ Σ n , deg(σ i g i ) ≤ d . If f ∈ Q(g ), we call the polynomials [σ 0 , σ 1 , . . . , σ k ] in an expression f = k i=0 σ i g i a SOS-certificate for f ∈ Q(g). The sets Q(g ) and Q(g )[d] depend on the polynomials g in the description of K. Denoted by Q(g) d := Q(g ) ∩ R[x] ≤d , by construction one has Q(g )[d] ⊂ Q(g ) d but, in
general, this inclusion is strict: in other words, for some g there exists a polynomial f ∈ Q(g ), of degree d, such that in any certificate f = k i=0 σ i g i , at least one product σ i g i has degree > d.

It is not even true that Q(g ) d ⊂ Q(g )[D] for possibly large D: these quadratic modules are called stable.

Definition 1 ([21, Sec. 4.1]) For d ∈ N, a quadratic module Q(g) is called stable in degree d if there exists D such that Q(g) d ⊂ Q(g)[D]. It is called stable if it is stable in every d ∈ N (we call the function D = D(d) a stability function for Q(g)).
Stability in a given degree depends on the generators g whereas stability depends only on the quadratic module Q(g ) and is equivalent to the existence of degree bounds for the representation f = i σ i g i ∈ Q(g ) that only depend on the degree of f , see for instance [START_REF] Netzer | Stability of quadratic modules[END_REF]. One example of stable quadratic module is the cone Σ n,2d = Q(0) 2d : indeed, every polynomial in Σ n,2d is a sum of squares of polynomials of degree at most d, that is,

Q(0) 2d = Q(0)[2d] and hence D(2d) = D(2d + 1) = 2d is a stability function for Σ n .
As well as Σ n,2d , truncated quadratic modules are semidefinite representable sets, that is linear images of feasible sets of SDP problems (also known as projected spectrahedra or spectrahedral shadows in the literature): given a description g for K = S(g ), computing one polynomial in Q(g )[D] amounts to solving a single SDP problem (cf. Section 3.3).

Definition 2 ([28])

A quadratic module Q(g) is called archimedean if there exists u ∈ Q(g) such that S(u) is compact.
Remark that if Q(g ) is archimedean, then S(g ) ⊂ S(u) thus S(g ) is compact. Archimedeanity and stability are often mutually exclusive properties, indeed, for n ≥ 2, an archimedean quadratic module is not stable.

Theorem 1 (Putinar's Positivstellensatz [START_REF] Putinar | Positive polynomials on compact semi-algebraic sets[END_REF]) Let K = S(g) be non empty, and assume Q(g) is archimedean. Then every polynomial positive on K belongs to Q(g).

The problem of bounding the degree D of the summands in a Putinar certificate f = k i=0 σ i g i for a polynomial f ∈ P(K) d , is called the effective Putinar Positivstellensatz, see [START_REF] Nie | On the complexity of Putinar's Positivstellensatz[END_REF][START_REF] Baldi | On the effective Putinar's Positivstellensatz and moment approximation[END_REF]. The work [START_REF] Baldi | On the effective Putinar's Positivstellensatz and moment approximation[END_REF] gives a bound for D as a function of d, n, of the polynomial f and of geometrical parameters of K, see [1, Th. 1.7], cf. Section 4.

Moments

A nonnegative Borel measure (below measure, for short) µ is a bounded nonnegative linear functional on the σ-algebra of Borel sets [START_REF] Basu | On the combinatorial and algebraic complexity of quantifier elimination[END_REF]. We say that µ satisfying (2) is a representing measure for the sequence y = (y α ) α∈N n d . If K is compact, the Stone-Weierstrass Theorem implies that a measure is uniquely determined by its (infinite-dimensional) sequence of monomial moments, see [START_REF] Marshall | Positive polynomials and sums of squares[END_REF]Cor. 3.3.1].

B(R n ). The support of µ is the complement of the largest open Borel set A ∈ B(R n ) such that µ(A) = 0, denoted by supp(µ). Let K ⊂ R n be a Euclidean closed set. A measure µ is supported on K if supp(µ) ⊂ K (in particolar it satisfies µ(R n \ K) = 0). For u ∈ K, we denote by δ u the Dirac measure supported on the singleton {u}. A finite linear combinations s i=1 c i δ u i of s Dirac measures is called s-atomic: its support is supp( s i=1 c i δ u i ) = {u 1 , . . . , u s }. For α ∈ N n , the (monomial) moment of exponent α of µ is the real number K x α dµ(x) as in
In this work, we address the following inverse problem for semialgebraic sets: Problem 1 (Truncated Moment Problem [START_REF] Schmüdgen | The moment problem[END_REF]) Let K ⊂ R n be a basic closed semialgebraic set. Given a finite sequence y = (y α ) α∈N n d of real numbers, determine whether y admits a representing measure supported on K.

For y = (y α ) α∈N n , the matrix M d (y) = (y α+β ) α,β∈N n d is called the moment matrix of order d of y. We recall that we denote by L y : R[x] ≤d → R the Riesz functional associated with y, defined by L y (x α ) = y α and extended linearly on R[x] ≤d .

Let now K = S(g ) be a basic closed semialgebraic set, and let

M (K) d = y = (y α ) α∈N n d ∈ R m : ∃ µ, supp(µ) ⊂ K, ∀α ∈ N n d , y α = K x α dµ(x)
denote the set of moments of order up to d of nonnegative Borel measures with support in K. The set M (K) d is in general not closed, as shown by the following example.

Example 1 ([7, Rem. 3.147]) Let n = 1, d = 4 and y = (1, 0, 0, 0, 1). Its second moment matrix M 2 (y) is positive semidefinite but the vector y is not representable by a nonnegative univariate measure, indeed y 2 = 0 but y 4 = 0. Thus y ∈ M (R) 4 . Nevertheless y ∈ M (R) 4 , the Euclidean closure of M (R) 4 : indeed y = lim →0 y where y = (1, 0, 2 , 0, 1) is the sequence of moments of degree ≤ 4 of the 3-atomic measure

µ = 4 2 δ1 + δ -1 + (1 -4 )δ 0
Let V ∨ be the dual vector space of a real vector space V , that is the set of R-linear functionals

L : V → R. If C ⊂ V is a convex cone, the set C * = {L ∈ V ∨ : L(a) ≥ 0, ∀a ∈ C} is called the dual cone of C. It is straightforward to see that M (K) d ⊂ P(K) * d
, and a non-trivial result is that equality holds for K compact:

Theorem 2 (Tchakaloff 's Theorem [START_REF] Tchakaloff | Formules de cubatures mécaniques à coefficients non négatifs[END_REF]) Let K ⊂ R n be compact and d ∈ N. Then

M (K) d = P(K) * d = {y ∈ R m : L y (p) ≥ 0, ∀p ∈ P(K) d }.
Theorem 2 is a finite-dimensional version of the Riesz-Haviland Theorem [START_REF] Curto | An analogue of the Riesz-Haviland theorem for the truncated moment problem[END_REF]. It is often used in the moment-SOS hierarchy, see [START_REF] Henrion | The moment-SOS hierarchy[END_REF]Lemma 1.7]. A modern statement and proof can be found e.g. in [19, Theorem 5.13], see also [START_REF] Blekherman | The core variety and representing measures in the truncated moment problem[END_REF].

An algorithm for the moment problem

We describe an algorithm based on semidefinite programming that solves Problem 1 for compact basic semialgebraic sets.

To do that, we first give a characterization of the interior of cones of nonnegative polynomials on compact sets K ⊂ R n (Section 3.1). Next we interpret Problem 1 as a conic feasibility problem and prove that its solvability is related to the feasibility type of the program (Section 3.2). Finally we describe our algorithm in Section 3.3.

Interior of P(K) d

Let K ⊂ R n be non-empty. For f ∈ R[x], we denote by f * := inf x∈K f (x), possibly -∞.

It is straightforward to construct examples of positive polynomials f ∈ P(K), on a noncompact set K, such that f ∈ Int(P(K)) even if f * > 0 (for instance 1 = lim →0 1 -x so 1 ∈ Int(P(R) 2 )). More generally, positive sum-of-squares polynomials of degree < d lie in the boundary of the cone of nonnegative polynomials of degree ≤ d (cf. [7, §4.4.3]).

The next folklore lemma shows that, for compact sets, the interior of P(K) d is exactly the set of positive polynomials over K.

Lemma 1 Let K ⊂ R n be non-empty, and let

d ∈ N. Then Int(P(K) d ) ⊂ {f ∈ R[x] ≤d : f * > 0}.
If K is compact, equality holds, and Int(P(K) d ) consists of exactly those polynomials in R[x] ≤d that are positive on K.

Proof: If f ∈ R[x] ≤d is such that f * ≤ 0, then (f -) * = f * -< 0 for all > 0, thus f -∈ P(K) d for all > 0, hence f ∈ Int(P(K) d ), which proves the sought inclusion ⊂. Now assume K is compact, and let f ∈ Int(P(K) d ). Then f is in the closure of the complement of P(K) d in R[x] ≤d , that is, f is the pointwise limit f = lim k→∞ f k of poly- nomials f k ∈ P(K) d , in particular, satisfying f * k < 0 for all k. Since K is compact, f * k = min x∈K f k (x) = f k (x k ) for some x k ∈ K. Let x ∈ K be a limit point of {x k } k ,
which exists by Bolzano-Weierstrass Theorem. Thus up to extracting a subsequence, one has 0

≥ lim k f * k = lim k f k (x k ) = f (x) ≥ f * , which shows the inclusion ⊃, thus the equality Int(P(K) d ) = {f ∈ R[x] ≤d : f * > 0}.
Since the infimum of a polynomial function on a compact set is its minimum, one has f * > 0 if and only if min x∈K f (x) > 0 if and only if f is positive on K, as claimed.

Remark 1 Theorem 1 and Lemma 1 ensure that if Q(g) is archimedean, then Int(P(K) d ) ⊂ Q(g) ∩ R[x] ≤d ⊂ P(K) d . Thus under this assumption, [1, Th. 1.7] yields a degree bound D = D(d, n, f, K) such that if f ∈ Int(P(K) d ) then f ∈ Q D (g). Since Q D (g) is semidefinite
representable, it can be sampled through semidefinite programming: solving such optimization problem yields an element of the boundary of Q D (g), thus this might not be sufficient to compute an element of Int(P(K) d ).

The following Corollary shows that one can get elements of Int(P(K) d ) as well through semidefinite programming, from the knowledge of a polynomial description g of K.

Corollary 1 Let g be such that Q(g) is archimedean, and let K = S(g). Let f ∈ Int(P(K) d ) and 0 < δ < f * = min K f . Then

1 δ f -1 ∈ Int(P(K) d ) and there exist σ δ 0 , σ δ 1 , . . . , σ δ k ∈ Σ n such that 1 δ f -1 = σ δ 0 + i σ δ i g i . Proof: The polynomial f -δ is positive on K, thus by Lemma 1, f -δ ∈ Int(P(K) d ) and hence (f -δ)/δ = 1 δ f -1 ∈ Int(P(K) d ), since Int(P(K) d ) is a cone. Since 1 δ f -1 is positive on K, we conclude by Theorem 1.
Corollary 1 can be rephrased as follows: if Q(g ) is archimedean and 0 < δ < f * = min K f , then f /δ ∈ 1 + Q(g ). Remark that (unless Q(g ) is stable in the degree of f ) the degrees of the SOS-multipliers for a SOS-certificate f /δ ∈ 1 + Q(g ) depend on δ and might be larger than the degrees for a SOS-certificate f ∈ Q(g ).

Moment problem as conic feasibility

Let C ⊂ V be a convex cone with non-empty interior, and let L ⊂ V be an affine space. The conic program associated with C and L is called Moreover, the following are equivalent and are strong alternatives to A 1 -A 2 -A 3 :

feasible if L ∩ C = ∅, otherwise infeasible. It is called strongly feasible if L ∩ Int(C) = ∅, and weakly feasible if it is feasible but not strongly. If L is a linear space (that is if 0 ∈ L) then {0} ⊂ (L ∩ C), thus L ∩ C is always feasible. If L is a hyperplane,
B 1 . y ∈ M (K) d ;
B 2 . The conic program L y ∩ P(K) d is strongly feasible.

Proof: The fact that A 1 and B 1 are strong alternatives is obvious, and the fact that the program L y ∩ P(K) d is always feasible (indeed, L y is linear) implies that A 2 and B 2 are strong alternatives. Hence we only have to prove the equivalence of A 1 , A 2 and A 3 .

We first prove that A 1 is equivalent to A 3 . For u ∈ K, denote by λ u = (u α ) α∈N n d ∈ M (K) d the sequence of moments of order ≤ d of the Dirac measure δ u . By [33, Th. 17.2], y admits a representing measure supported on K, if and only if it admits a representing atomic measure

µ = s i=1 c i δ u i , where s ≤ dim R[x] ≤d = n+d d
and for some c i > 0 and u 1 , . . . , u s ∈ K. For every p ∈ R[x] ≤d , one deduces

L y (p) = K p dµ = s i=1 c i p(u i )
and thus p ∈ L y if and only if s i=1 c i p(u i ) = 0: then for p ∈ P(K) d , we conclude that p must vanish on {u 1 , . . . , u s }. We deduce that A 1 and A 3 are equivalent.

We prove now that A 1 and A 2 are equivalent. By Theorem 2, we know that A 1 holds if and only if L y is non-negative over P(K) d : this is the case if and only if the cone P(K) d is contained in the closed half-space L + y = {p ∈ R[x] ≤d : L y (p) ≥ 0}, and since 0 ∈ L y , this is equivalent to L y being a supporting hyperplane and the program being weakly feasible.

When the conic program in Proposition 1 is weakly feasible, the set L y ∩ P(K) d is an exposed and proper face of P(K) d defined by vanishing on the finite set defined in Item A 3 . See also [START_REF] Blekherman | Dimensional differences between faces of the cones of nonnegative polynomials and sums of squares[END_REF] and [START_REF] Blekherman | Semidefinite optimization and convex algebraic geometry[END_REF]Sec. 4.4]. On the contrary, if the conic program is strongly feasible, we give the following definition.

Definition 3 Let y ∈ M (K) d . A polynomial p ∈ L y ∩ Int(P(K) d ) is called a unrepresen- tability certificate for y in K.
Corollary 2 shows how to compute explicit unrepresentability certificates for Problem 1.

Corollary 2 Assume that Q(g) is archimedean and that conditions B 1 -B 2 of Proposition 1 hold. There exists p ∈ Int(P(K) d ) such that p ∈ 1 + Q(g), L y (p) = 0 and p * > 1 is arbitrarily large.

Proof: Property B 2 of Proposition 1 ensures that there exists f ∈ L y ∩ Int(P(K) d ), that is f is positive over K and L y (f ) = 0. Let 0 < δ < f * , and let p = 1 δ f ∈ Int(P(K) d ). From Corollary 1 we get that p -1 ∈ Q(g ), that is, p ∈ 1 + Q(g ). Moreover L y (p) = 1 δ L y (f ) = 0, and p * = f * /δ > 1 is arbitrarily large.

Remark 2 If y ∈ Q m , then p in Corollary 2 can be chosen with rational coefficients. Indeed, the hyperplane L y is defined by an equation with rational coefficients, and L y ∩ Int(P(K) d ) is a non-empty open subset of L y , hence it contains a rational point f . Choosing δ ∈ Q in the proof of Corollary 2 is thus sufficient to get a rational certificate.

Nevertheless, let us recall from [START_REF] Scheiderer | Sums of squares of polynomials with rational coefficients[END_REF] that there exist polynomials in Q[x], that are sums of squares as elements of R[x] but not as elements of Q[x]. In our context, this means that the rational unrepresentability certificate p might not admit rational certificates of positivity showing that p ∈ 1 + Q(g) (see also [START_REF] Naldi | Conic programming: infeasibility certificates and projective geometry[END_REF] for the existence of rational certificates in conic programming).

Any polynomial p as in Corollary 2 is such that p -1 ∈ Q(g ). In particular, there exists The following example shows that the certificate of Corollary 2 might exist in non-archimedean contexts.

D > 0 such that p -1 ∈ Q(g )[D], that is p ∈ 1 + Q(g

Example 2

The vector y = (1, 1, 0) is not a univariate moment vector (indeed the moment matrix M 1 (y) = ( 1 1 1 0 ) is not positive semidefinite). Remark that R = S(0) and Q(0) = Σ n is not archimedean, but an unrepresentability certificate in the spirit of Corollary 2 exists. Indeed L y (p 0 + p 1 x + p 2 x 2 ) = p 0 + p 1 and the following identity holds for all p ∈ L y :

p 0 -p 0 x + p 2 x 2 = 1 + 1 x p 0 -1 -p 0 2 -p 0 2 p 2 1 x .
The identity is a SOS-certificate for p ∈ 1 + Q(g) if and only if the Gram matrix on the right hand side is positive semidefinite. This yields a spectrahedral representation of the set of unrepresentability certificates of y:

p ∈ R[x] ≤2 : p 1 = -p 0 , p 0 -1 - p 0 2 - p 0 2 p 2 0 For instance the polynomial p = 2 -2x + x 2 = 1 + (1 -x) 2 ∈ (1 + Q(0)) ∩ L y .
We remark that the existence of such p implies y ∈ M (R) 2 . Indeed, if y ∈ M (R) 2 , then by A 3 of Proposition 1, there would exist a representing measure with finite support, thus one would have The polynomial p in Example 2, together with its positivity certificate 1 + (1 -x) 2 , allows to check rigourously the unrepresentability of y = (1, 1, 0). Nevertheless, in general the archimedianity hypothesis cannot be dropped, as shown in Example 3.

y ∈ M ([-R, R]) 2 for some R > 0. Now, since 1 + Q(0) ⊂ 1 + Q(R 2 -x 2 ) for every R > 0, we deduce that p ∈ 1 + Q(R 2 -x 2 ),
Example 3 Let y = (1, 0, 0, 0, 1) be the vector of Example 1, with g = 0. The semidefinite program in Corollary 2 is infeasible: indeed Q(0) = Σ n and it is easy to check that there is no polynomial p = 4 i=0 p i x i satisfying L y (p) = p 0 + p 4 = 0

p = 1 + 1 x x 2 X 1 x x 2 T X 0.
Indeed, the constraints imply that p 4 = -p 0 = -(1 + X 11 ) < 0, thus p is negative at infinity, in particular, p ∈ P(R).

Nevertheless for every R > 0, with g = (R -x, R + x), Corollary 2 ensures that there exists

p R ∈ 1 + Q(R -x, R + x) such that L y (p R ) = 0, that certifies that y ∈ M ([-R, R]) 4 :
the polynomial p R is any solution of the following parametric linear matrix inequality with σ 0 , σ 1 , σ 2 ∈ Σ 1 :

p 0 + p 1 x + p 2 x 2 + p 3 x 3 -p 0 x 4 = 1 + σ 0 + σ 1 (R -x) + σ 2 (R + x).
The fact that g = (R -x, R + x) is the natural description of [-R, R] (see [START_REF] Marshall | Positive polynomials and sums of squares[END_REF]Sec. 2.7]) implies that Q(g) is stable, with stability function D(d) = d (see e.g. [START_REF] Schmüdgen | The moment problem[END_REF]Prop. 3.3]), and hence one can assume σ 0 ∈ Σ 1,4 and σ 1 , σ 2 ∈ Σ 1,2 .

We terminate this series of examples with a bivariate one. We claim that there is no nonnegative Borel measure supported on the unit ball

K = {a = (a 1 , a 2 ) ∈ R 2 : a 2 1 + a 2 2 ≤
1} whose moments up to degree 6 agree with y. Remark that the semialgebraic set

K = S(1 -x 2 1 -x 2 2 ) is compact and Q(1 -x 2 1 -x 2 2
) is archimedean, in particular it is not stable. We give below in Example 5 an unrepresentability certificate of small degree certifying that y ∈ M (K) d , proving our claim.

SDP-based algorithm

The results of Section 3.2 yield the following alternatives for Problem 1. Given y = (y α ) α∈N n d and K = S(g ), according to Proposition 1: either y ∈ M (K) d , in which case a certificate of unrepresentability is given by a polynomial p ∈ Int(P(K) d ) such that L y (p) = 0 (Corollary 2) or y ∈ M (K) d , in which case there exists an atomic measure µ = s i=1 c i δ u i representing y (Property A 3 in Proposition 1). We describe our main algorithm.

Algorithm 1 certify moment

Input: n, d ∈ N A vector y ∈ R m , with m = n+d d Polynomials g = (g 1 , . . . , g k ) ∈ R[x] k A threshold D ∈ N Output:
Either (p, Σ) where p ∈ R[x] ≤d satisfies Corollary 2, and

Σ ∈ Σ k+1 n is a certificate for p ∈ 1 + Q(g)[D]
or a measure µ = c i δ u i satisfying A 3 in Proposition 1 1: (p, Σ) ← find certificate(n, d, y, g, D) 2: if p = [] then return (p, Σ) 3: return find measure(n, d, y, g) Algorithm 1 depends on two subroutines. The first one, find certificate at Step 1, returns, if it exists, an unrepresentability certificate p ∈ R[x] ≤d for y, together with a SOS-certificate for p as element of 1 + Q(g )[D].

Algorithm 2 find certificate

1: p ← [], Σ ← [], g 0 ← 1 2: Find p ∈ R[x] ≤d and σ 0 , σ 1 , . . . , σ k ∈ Σ n such that L y (p) = 0 p = 1 + k i=0 σ i g i , deg(σ i g i ) ≤ D 3: Σ ← [σ 0 , σ 1 , . . . , σ k ] 4: return (p, Σ)
Algorithm 2 can be performed by solving one (finite-dimensional) SDP feasibility program whose unknowns are p, σ 0 , σ 1 , . . . , σ k . First, the constraint L y (p) = 0 is linear in p. Next, denoting

δ i = (D -deg(g i ))/2 , the constraints σ i ∈ Σ n and deg(σ i g i ) ≤ D are equivalent to the existence of a symmetric matrix X i 0, of size δ i +n n , such that σ i = v T X i v, where v is a linear basis of R[x] δ i . Finally the constraint p = 1 + Q(g )[D] in
Step 2 is affine linear in p and in the entries of X 0 , X 1 , . . . , X k . We give upper bounds for the value of D in Section 4.

The second routine, find measure, is called if and only if Algorithm 1 reaches Step 3. It returns a s-atomic measure representing the vector y, for some s ≤ n+d d .

Algorithm 3 find measure

1: s = 1 2: while s ≤ n+d d do 3: Find c ∈ R s and U = (u 1 , . . . , u s ) ∈ (R n ) s s.t. u 1 , . . . , u s ∈ S(g) c 1 u α 1 + • • • + c s u α s = y α for all α ∈ N n d 4:
if solution exists then return (c, U )

s ← s + 1

The routine find measure can be performed by existing algorithms computing one point per connected component of basic semialgebraic sets applied to the set of elements (c, U ) ∈ R s ×(R n ) s satisfying the inequalities defining K and the polynomial equations

i c i u α i = y α , α ∈ N n d .
The equations have a multivariate Vandermonde structure. A precise complexity analysis of find measure is left to future work.

We define now an integer function of the input of Problem 1. Definition 4 Let K = S(g) and y ∈ M (K) d . The unrepresentability degree of y in K is the minimum integer D = D(n, d, y, g) such that there exists

p ∈ 1 + Q(g)[D] satisfying L y (p) = 0.
The unrepresentability degree is well defined, according to Corollary 2. We prove the correctness of Algorithm 1.

Theorem 3 (Correctness) Let y = (y α ) α∈N n d , and let g = (g 1 , . . . , g k ) ∈ R[x] k be such that Q(g) is archimedean. There exists D = D(n, d, y, g) ∈ N such that Algorithm 1 with input (n, d, y, g, D) terminates and is correct.

Proof: Let K = S(g ). Since Q(g ) is archimedean, K is compact and by Lemma 1, Int(P(K) d ) = {p ∈ R[x] ≤d : p(a) > 0, ∀ a ∈ K}. By Theorem 1, if a polynomial is positive on K, then it belongs to Q(g ). We deduce that Int(P(K) d ) ⊂ Q(g ) ∩ R[x] ≤d ⊂ P(K) d .
We claim that for D ∈ N large enough, then y ∈ M (K) d if and only if certify moment returns a polynomial p at Step 2, that is, if and only if the semidefinite program at Step 1 is feasible. For the reverse implication, suppose that the semidefinite program at Step 1 is feasible for some degree D. Let p ∈ R[x] ≤d be a solution of such program. Then p = 1 + q for some q ∈ Q(g )[D], in particular, p * ≥ 1 on K, that is, p is positive on K. Since L y (p) = 0, and by compactness of K, one has p ∈ Int(P(K)) ∩ L y and again applying Proposition 1 one concludes y ∈ M (K) d . This proves the claim.

Finally, remark that this also shows that y ∈ M (K) d if and only if certify moment reaches Step 3. If this is the case, find measure computes the support and the weights of an s-atomic measure representing y: such measure exists for some s ≤ n+d d , according to Proposition 1.

Bound on the unrepresentability degree

A priori bounds on the degree of Putinar certificates that only depend on the degree of the polynomial exist for stable quadratic modules. Nevertheless, stability and archimedeanity properties are mutually exclusive in dimension n ≥ 2.

In this section we give general bounds for the unrepresentability degree of a vector y ∈ M (K) d in a compact basic semialgebraic set K = S(g ). A key ingredient to do this is the use of already existing quantitative analysis of computer algebra algorithms performing quantifier elimination over the reals.

We first recall the following bound on the degree of a Putinar representation of a polynomial in Int(P(K) d ), for an archimedean quadratic module Q(g ), given in [START_REF] Baldi | On the effective Putinar's Positivstellensatz and moment approximation[END_REF]. Let Under the following assumptions:

1 -i x 2 i ∈ Q(g ) g i ≤ 1
2 for all i = 1, . . Let now (n, d, y, g ) be the input of Algorithm 1. We assume in the whole section that y ∈ Q N , with N = n+d d , that the basic semialgebraic set K = S(g ) ⊂ R n is defined by polynomial inequalities g = (g 1 , . . . , g k ) ⊂ Q[x] k and that Q(g ) is archimedean. In the following, we provide an upper bound B on p : this will yield a lower bound on (p -1), hence an upper bound on D.

As already said, in order to do that, we consider the formulation in terms of quantifier elimination over the reals that solves the truncated moment problem in the sense of Proposition 1 and Corollary 2, and we use quantitative results on quantifier elimination over the reals from [START_REF] Basu | Algorithms in real algebraic geometry[END_REF].

Consider the following formula with quantified variables x = (x 1 , . . . , x n ) and parameters the unknown coefficients

{p α : α ∈ N n d } of a polynomial p ∈ R N : L y (p) = 0 ∧ ∀x ∈ R n k i=1 g i (x) ≥ 0 ⇒ p(x) > 1 . (4)
Lemma 2 Let S ⊂ R N be the semialgebraic set defined by the quantifier-free formula obtained from (4) after eliminating the quantified variables. Then S is an open subset of L y with respect to the induced topology.

Proof: For p ∈ R[x] ≤d and A ⊂ R[x] ≤d , we denote by p + A := {p + q : q ∈ A}. Since Q(g ) is archimedean, K is compact and according to Lemma

= (-1 + L y ) ∩ Int(P(K) d ) is open in -1 + L y and hence S = 1 + T is open in 1 + (-1 + L y ) = L y .
Observe that the constraint L y (p) = 0 in ( 4) is linear in p and y = 0. Thus one of the coefficients of p, say p α , can be eliminated, yielding an formulation of the quantifier elimination problem (4):

∀x ∈ R n k i=1 g i (x) ≥ 0 ⇒ p(x) > 1 ( 5
)
where p is the polynomial obtained when substituting p α in p by a linear form in the other coefficients p α using L y (p) = 0. Below we abuse of notation and consider the set S defined in Lemma 2 and by the previous formulae as embedded in L y identified with R N -1 . We conclude that the set S has non-empty interior in R N -1 .

We denote by

d g = max i deg(g i ), i = 1, . . . , k
and by τ g the maximum bit size of the coefficients of the g i 's. Note that we can multiply in (5) the polynomials g i by the (positive) least commun multiple of the denominators of their coefficients to obtain equivalent inequalities but with coefficients in Z. These least common multiples have height bounded by ( n+dg n + 1)τ g .

Further, we denote by τ y the maximum bit size of the coefficients of y. As above, the equation L y (p) = 0 can be rewritten with coefficients in Z of bit size bounded by (N + 1)τ y . We set

τ = max (N + 1)τ y , n + d g n τ g . (6)
Note that τ is a bound on the integer coefficients of the polynomial constraints in (5) once we have multiplied each of them by the least common multiple of the denominators of their coefficients.

Finally, let δ = max(d+1, d g ). Note that δ dominates the maximum degree of the polynomial constraints in (5), indeed, the polynomial p ∈ R[x, p α : α ∈ N n d ] has degree d in x and has degree 1 with respect to its unknown coefficients. Proof: We start by providing some quantitative bounds on the quantifier-free formula which defines S ⊂ R N -1 obtained by eliminating the quantified variables x = (x 1 , . . . , x n ) in [START_REF] Blekherman | Dimensional differences between faces of the cones of nonnegative polynomials and sums of squares[END_REF]. By [START_REF] Basu | Algorithms in real algebraic geometry[END_REF]Theorem 14.16] such a formula satisfies the following properties It can be obtained with polynomials of degree lying in (δ + 1) O(n) ; the bit size of the coefficients of these polynomials lies in τ (δ + 1) O(n(N -1)) ; this formula is a disjunction of

k n+1 δ O(n) conjunctions of k n+1 δ O(n) disjunctive formulas of polynomial inequalities involving k n+1 δ O(n) polynomials.
Recall that since S is open, it coincides with its interior Int(S). We aim at computing one point with rational coordinates in S. To do this, we just put these disjunctions in closed form, replace non-strict inequalities by strict inequalities and call an algorithm for computing at least one point with rational coordinates in S; see e.g. [START_REF] Basu | On the combinatorial and algebraic complexity of quantifier elimination[END_REF]Theorem 4.1.2].

Note that the input to such an algorithm is a system of polynomial strict inequalities in R[p α : α ∈ N n d ] of degree (δ + 1) O(n) with bit size coefficients in τ (δ + 1) O(n(N -1)) . By [2, Theorem 4.1.2], if the semialgebraic set defined by the input is non-empty, then it outputs a point with rational coordinates of bit size bounded by τ O(1) (k(δ + 1)) O(n(N -1)) .

All in all, this bounds the bit size of the coefficients of some polynomial p ∈ S (with rational coordinates). Since the number of these coefficients is N -1, the 2-norm of p still lies in τ O(1) (k(δ + 1)) O(n(N -1)) .

Finally, observe that p = min a∈[-1,1] n p(a) is bounded above by the 2-norm of p.

Corollary 3 Let τ y and τ g bound the bit-size of y and g, respectively, and let d g be a bound on the degrees of g 1 , . . . , g k . Let τ be as in [START_REF] Blekherman | The truncated K-moment problem for closure of open sets[END_REF] We terminate with a bivariate example showing that the bound of Corollary 3 is usually quite pessimistic. ) is compact, p certifies that the conic program defined in Proposition 1 is strongly feasible, in other words, that y ∈ M (K) d .

Conclusions and perspectives

The goal of this work was to undertake a systematic analysis of the computational complexity of the truncated moment problem on semialgebraic sets. Preliminary results concern the existence of algebraic certificates for vectors that are not representable as moments of measures, and upper bounds on the degree of SOS representations of these certificates.

Our contribution offers several challenges and research directions in the computational aspects of the truncated moment problem, let us mention a few. One of these is the need of efficient algorithms for classes of polynomial systems with Vandermonde structure as that defined in the routine find measure.

A second one is to refine the quantifier-elimination bound given in Corollary 3. Unlike the viewpoint of the so-called effective Putinar Positivstellensatz introduced in [START_REF] Baldi | On the effective Putinar's Positivstellensatz and moment approximation[END_REF], for which degree bounds depend on the polynomial itself, it is clear from our analysis that for the complexity analysis of the TMP one needs to give uniform degree bounds that only depend on the input of the TMP. One way of getting such uniform bounds is to consider manifestly positive polynomials such as those in 1 + Q(g ) for compact K = S(g ).

A final perspective is to extend our analysis to the more general case of basic closed semialgebraic sets, not necessarily compact (see e.g. [START_REF] Blekherman | The truncated K-moment problem for closure of open sets[END_REF] and [4, ]).

Proposition 1

 1 then the corresponding program is weakly feasible if and only if L ∩ C is a proper face of C and in this case L is called a supporting hyperplane for C: geometrically, C is contained in one of the two closed half-spaces bounded by L, and L is tangent to the boundary of C. Let y = (y α ) α∈N n d ∈ R m , m = n+d d with y 0 > 0. Let L y ∈ (R[x] ≤d ) ∨ be the Riesz functional of y, and L y = {p ∈ R[x] ≤d : L y (p) = 0}. Let K = S(g) ⊂ R n . The following are equivalent:A 1 . y ∈ M (K) d ; A 2 . The conic program L y ∩ P(K) d isweakly feasible;A 3 . There exist u 1 , . . . , u s ∈ K, with s ≤ m, such that y admits a representing measure µ with supp(µ) = {u 1 , . . . , u s } and L y ∩ P(K) d = {p ∈ P(K) d : p(u 1 ) = 0, . . . , p(u s ) = 0}.

  )[D]. Bounds for D are given, e.g., in [1, Th. 1.7]. Below we give a bound on D as a function of the input of Algorithm 1, see Definition 4 and Section 4.

  thus p ∈ Int(P([-R, R]) 2 ) (according to Lemma 1) and thus y ∈ M ([-R, R]) 2 , for every R > 0 (by B 2 , Proposition 1).

Example 4

 4 Let n = 2, d = 6 and let y = (y α ) α∈N 2 6 be the vector in R 28 whose non-zero entries are y 00 = 32 y 22 = 30 y 20 = y 02 = 34 y 60 = y 06 = 128 y 40 = y 04 = 43 y 42 = y 24 = 28.

  Assume y ∈ M (K) d and let D the unrepresentability degree of y in K. By Corollary 2, there exists a polynomial p ∈ 1 + Q(g )[D] ⊂ Int(P(K) d ), such that L y (p) = 0. In other words, find certificate returns (p, Σ) with p = [], and hence Algorithm 1 returns its output at Step 2.

f

  ∈ Int(P(K) d ). We denote by (f ) := f * / f where f * = min a∈K f (a) and f = max a∈[-1,1] n f (a).

  . , k then by [1, Th. 1.7] there exists a function γ = γ(n, g ) such that f ∈ Q(g )[D] for D of the order of γ(n, g ) d 3.5 Ln (f ) -2.5 Ln (3) where c and L are the Lojasiewicz coefficients as they are defined in [1, Def. 2.4].

  Assume y ∈ M (K) d , and let p ∈ 1 + Q(g ) satisfy L y (p) = 0 and p * = 1 + m for some arbitrary constant m > 0. Such unrepresentability certificate exists by Corollary 2. From (3), one has that p ∈ 1 + Q(g )[D] with D depending on (p -1) = p * -1 p -1 ≥ 1 1 + p where the last inequality derives from Corollary 2 choosing without loss of generality m = 1.

Proposition 2

 2 There exists p ∈ S with p ≤ B for B in τ O(1) (k(δ + 1)) O(n(N -1)) .

  and δ = max{d + 1, d g }. If y ∈ M (K) d , then the degree of unrepresentability of y in K is in γ(n, g) d 3.5 Ln τ O( Ln) (k(δ + 1)) O( Ln 2 (N -1)) . Proof: With the notation introduced in (3), by Corollary 2, there exists p ∈ 1 + Q(g )[D] and by applying [1, Th. 1.7] and Proposition 2, the degree D is bounded above by γ(n, g ) d 3.5 Ln (p -1) -2.5 Ln ≤ γ(n, g ) d 3.5 Ln (1 + p ) 2.5 Ln ≤ γ(n, g ) d 3.5 Ln (τ O(1) (k(δ + 1)) O(n(N -1)) ) 2.5 Ln = γ(n, g ) d 3.5 Ln τ O( Ln) (k(δ + 1)) O( Ln 2 (N -1))

Example 5 ( 9 ( 1 -x 2 1 -x 2 2 )

 59122 Example 4 continued) Let y ∈ R 28 be the vector defined in Example 4. Consider the polynomialp = 1 + 8 One checks that p ∈ 1 + Q(1 -x 2 1 -x22 ) and L y (p) = y 00 (1 + 8 9 ) -8 9 (y 20 + y 02 ) = 0. Since K = S(1 -x 2 1 -x 2 2

  1, Int(P(K) d ) consists exactly of polynomials that are positive on K. Thus S = L y ∩ (1 + Int(P(K) d )). Moreover by Proposition 1, we know that L y ∩ Int(P(K) d ) is non-empty and open in L y . Now since L y is a linear space intersecting Int(P(K) d ), which is an open convex cone, so does the affine space -1 + L y . Thus T

Acknowledgments

The authors thank Lorenzo Baldi for discussions concerning the degree bounds of SOScertificates. This work is supported by the European Commission Marie Sklodowska-Curie Innovative Training Network POEMA (Polynomial Optimization, Efficiency through Moments and Algebra, 2019-2023); by the Agence Nationale de la Recherche (ANR), grant agreements ANR-18-CE33-0011 (SESAME), ANR-19-CE40-0018 (De Rerum Natura), ANR-21-CE48-0006-01 (HYPERSPACE); by the joint ANR-Austrian Science Fund FWF grant agreement ANR-19-CE48-0015 (ECARP); by the EOARD-AFOSR grant agreement FA8665-20-1-7029.