Edouard Pauwels

Didier Henrion

Jean-Bernard Lasserre

LINEAR CONIC OPTIMIZATION FOR INVERSE OPTIMAL CONTROL *

Keywords: f, X, U Unknown: µ, v Output: Optimal measure µ Inverse control Controled trajectories viewed as occupation measure µ Conic duality: Data: µ Unknown: l, v Output: Lagrangian l

We address the inverse problem of Lagrangian identification based on trajectories in the context of nonlinear optimal control. We propose a general formulation of the inverse problem based on occupation measures and complementarity in linear programming. The use of occupation measures in this context offers several advantages from the theoretical, numerical and statistical points of view. We propose an approximation procedure for which strong theoretical guarantees are available. Finally, the relevance of the method is illustrated on academic examples.

1. Introduction. In the context of nonlinear optimal control, we are interested in the inverse problem of Lagrangian identification from given trajectories. This identification should be carried out such that solving the direct optimal control problem with the identified Lagrangian would allow to recover the given trajectories.

Inverse problems of calculus of variations are old topics that have attracted a renewal of interest in the context of optimal control, especially in humanoid robotics [START_REF] Arechavaleta | An optimality principle governing human walking[END_REF]. Relevant aspects of the problem are not well understood and many issues still need to be addressed to propose a tool that could be used in experimental settings. The work presented here constitutes a step in this direction. A preliminary conference version [START_REF] Pauwels | Inverse optimal control with polynomial optimization[END_REF] originally introduced our optimization framework as a tool to solve the inverse problem numerically. The current paper extends this work in many ways. In particular, by using the (quite general) concept of occupation measures we can propose a broad definition of inverse optimality and we also rigorously justify most of the approximations behind the numerical results reported in [START_REF] Pauwels | Inverse optimal control with polynomial optimization[END_REF]. Many aspects of this work parallel the results of [START_REF] Lasserre | Nonlinear optimal control via occupation measures and LMI relaxations[END_REF] about direct optimal control with polynomial data.

1.1. Motivation. The principle of optimality (or stationarity) is very important as a conceptual tool to describe laws of phenomenon are observed in nature (e.g. Fermat's principle in optics, Lagrangian dynamics in mechanics). Beyond physics, similar tools and arguments are used to describe and model the behaviour of living systems in biology [START_REF] Rosen | Optimality principles in biology[END_REF] or decision making agents in economics [START_REF] Kamien | Dynamic optimization: the calculus fo variations and optimal control in economics and management[END_REF]. Of more important interest to us is the application of the optimality principle to model the motion of living organisms [START_REF] Todorov | Optimality principles in sensorimotor control[END_REF]. In our technological context, this constitutes a hot topic. Promising expectations for these types of model include:

• The conceptual understanding of general laws that govern decision taking processes related to living organism motion, including human motion [START_REF] Arechavaleta | An optimality principle governing human walking[END_REF]. • The ability to use these general laws to reproduce and synthetise motion behaviours for new tasks with unknown space configuration. In this context, the principle of optimality only constitutes one possible conceptual tool to understand motion. There is a debate regarding its validity [START_REF] Friston | What is optimal about motor control?[END_REF] or its direct applicability in robotics applications [START_REF] Laumond | Optimality in robot motion: optimal versus optimized motion[END_REF]. These illustrate the fact that this idea constitutes an active subject of research, with a strong connextion with applications. Direct optimal and inverse optimal control flow chart. System description is given by the dynamics f , the state constraint set X and control constraint set U which are all assumed to be fixed. We emphasize that the Lagrangian and the occupation measure have symmetric roles for the direct and inverse problems. We also note that the output of the inverse problem is a Lagrangian. Solving the direct optimal control problem for new initial conditions is an important question but remains secondary regarding inverse optimality which is the focus of this work.

In many situations however, the cost related to the motion of a system is unknown or does not correspond to direct intuition. In these cases, as clearly emphasized in [START_REF] Todorov | Optimality principles in sensorimotor control[END_REF]: "It would be very useful to have a general data analysis procedure that infers the cost function given experimental data and a biomechanical model". Our contribution is to investigate the mathematical meaning of "inferring cost function from data" and we propose a numerical method to address problems of this type based on inverse optimality. We emphasize that this paper is "only" concerned with this question. In particular we do not address the issue of interpreting the inferred cost function or solving direct problems for new unseen conditions. We solely focus on the task of inferring a cost function from data. This constitutes a nontrivial shift in term of point of view compared to usual questions arising when dealing with direct optimal control problems. We hope to convince the reader that there are crucial differences between inverse and direct optimal control and that it is worth investigating the former within an appropriate context with somewhat different questions in mind.

The backbone of the proposed approach and its relation with the direct problem of optimal control is presented in Figure 1.1. It is important to understand the symmetric role of the Lagrangian and the occupation measure representing the input trajectories. As a matter of fact, since the input of the inverse problem is a set of trajectories (supposedly optimal for a certain Lagrangian), many aspects of the existence of minimizers that are crucial in direct optimal control, are not relevant for inverse problems since the "optimal" trajectories are given. For example, there is no need to recompute optimal trajectories for direct problems with initial conditions already considered in the input data since by inverse optimality, the input trajectories are optimal with respect to the identified Lagrangian.

1.2. Context. Since its introduction by Kalman [START_REF] Kalman | When is a linear control system optimal[END_REF], the inverse problem of optimal control has been studied in linear settings [START_REF] Anderson | Linear optimal control[END_REF][START_REF] Jameson | Inverse problem of linear optimal control[END_REF][START_REF] Fujii | A complete optimality condition in the inverse problem of optimal control[END_REF][START_REF] Nori | Linear optimal control problems and quadratic cost functions estimation[END_REF] leading to many nonlinear variations [START_REF] Thau | On the inverse optimum control problem for a class of nonlinear autonomous systems[END_REF][START_REF] Moylan | Nonlinear regulator theory and an inverse optimal control problem[END_REF][START_REF] Casti | On the general inverse problem of optimal control theory[END_REF][START_REF] Freeman | Inverse optimality in robust stabilization[END_REF]. In these works the input data of the problem is a characteristic of a class of trajectories often given in the form of a control law. This contrasts with the setting we propose to study, for which the input is a set of trajectories which could come from physical experiments. This motivates the work of [START_REF] Chittaro | On inverse optimal control problems of human locomotion: stability and robustness of the minimizers[END_REF] and [START_REF] Ajami | How humans fly[END_REF] about well-posedness of the inverse problem, both in the context of unicycle dynamics in robotics and strictly convex positive Lagrangians.

On the other hand, to treat the inverse problem several authors have proposed numerical methods based on the ability to solve the direct problem [START_REF] Mombaur | From human to humanoid locomotion-an inverse optimal control approach[END_REF], also in the context of Markov decision process [START_REF] Abbeel | Apprenticeship learning via inverse reinforcement learning[END_REF][START_REF] Ratliff | Maximum margin planning[END_REF] or based on a discretized version of the direct problem [START_REF] Puydupin-Jamin | A convex approach to inverse optimal control and its application to modeling human locomotion[END_REF][START_REF] Keshavarz | Imputing a convex objective function[END_REF].

Our approach is different and based on occupation measures, an abstract and quite general tool to handle trajectories (and their weak limits) of feasible solutions of classical control problems. Formulating the (direct) control problem on appropriate spaces of measures amounts to relaxing the original problem. In most applications, both relaxed and original problems have same optimal value [START_REF] Vinter | Convex duality and nonlinear optimal control[END_REF][START_REF] Vinter | The equivalence of strong and weak formulations for certain problems in optimal control[END_REF][START_REF] Gaitsgory | Linear Programming Approach to Deterministic Infinite Horizon Optimal Control Problems with Discounting[END_REF]. However the relaxed formulation has the crucial advantage that compactness holds in a certain weak sense: As a matter of fact, many optimization problems over appropriate spaces of measures attain their optimum, whereas most optimization problems over smaller functional spaces (e.g. continuous functions, or Lebesgue integrable functions) typically have no optimal solution. At last but not least, for control problems with polynomial data, the relaxed problem can be formulated as an optimization problem on moments of occupation measures. By combining this with relatively recent advances in real algebraic geometry [START_REF] Putinar | Positive polynomials on compact semi-algebraic sets[END_REF] and in numerical optimization [START_REF] Lasserre | Moments, positive polynomials and their applications[END_REF] one may thus provide a systematic numerical scheme to approximate effectively relaxed solutions of optimal control problems [START_REF] Lasserre | Nonlinear optimal control via occupation measures and LMI relaxations[END_REF].

1.3. Contribution. We choose the setting of free terminal time optimal control which is consistent with many physical experiments that one can think of. But the same approach with ad hoc modifications is also valid in the fixed terminal time setting.

• In our opinion, occupation measures are the perfect abstract tool to formally express the fact that we consider a (possibly uncountably infinite) superposition of trajectories as input data of the inverse control problem. We then propose a general formulation of the inverse problem based on occupation measures and complementarity in linear programming. A relaxation of the well known Hamilton-Jacobi-Bellman (HJB) sufficient optimality condition appears in our formulation as for the usual direct optimal control problem [START_REF] Hernández-Hernández | The linear programming approach to deterministic optimal control problems[END_REF]. This formulation is shown to be consistent with what is commonly expected regarding inverse optimality.

It is worth noting that when using the HJB optimality conditions, the situation is completely symmetric for the direct and inverse control problems. In both cases the HJB optimality conditions are used to certify the global optimality of trajectories. But in the former the Lagrangian is known and HJB provide conditions on the optimal state-control trajectories (to be determined) whereas in the latter the "optimal" statecontrol trajectories are known and HJB provide conditions on the Lagrangian (to be determined) for the given trajectories to be optimal. (In both cases the optimal value function is considered as an auxiliary "variable".)

• Furthermore, this framework allows to further characterize the space of solutions associated with a given inverse optimal control problem. This viewpoint is different from what has been proposed in previous (theoretical and numerical) contributions to this problem [START_REF] Mombaur | From human to humanoid locomotion-an inverse optimal control approach[END_REF][START_REF] Puydupin-Jamin | A convex approach to inverse optimal control and its application to modeling human locomotion[END_REF][START_REF] Chittaro | On inverse optimal control problems of human locomotion: stability and robustness of the minimizers[END_REF][START_REF] Ajami | How humans fly[END_REF] which, implicitly or explicitly, involve strong (and, in our opinion, overly restrictive) constraints on the class of functions in which the candidate Lagrangians are searched.

• The weak formulation of direct optimal control problems via occupation mea-sures is elegant and powerful but also involves difficult technical questions regarding potential gaps between classical and generalized control problems. Using inverse optimality, we justify a posteriori that this discussion can be partially mitigated for the inverse problem. This striking difference between direct and inverse problems is due to the symmetric roles of the Lagrangian and occupation measure and the fact that the occupation measure is given and fixed for the inverse problem.

• Remarkably, despite the abstract setting of occupation measures, the proposed formulation is amenable to explicit numerical approximations via a hierarchy of semidefinite programs 1 . Indeed in the context of polynomial dynamics and semi-algebraic constraints, both the optimal value function and Lagrangian used in the (relaxed) HJB optimality conditions can be approximated with polynomials. We show that such a reinforcement is coherent in the sense that no polynomial solution to the inverse problem is lost.

• Finally, in usual experimental settings one does not have access to complete trajectories. Instead one is rather given finitely many data points sampled from trajectories. But results from probability applied to our occupation measures allow to formalize the fact that we only work with "samples". In addition, in this framework one may use empirical processes and statistical learning theory [START_REF] Vapnik | An overview of statistical learning theory[END_REF][START_REF] Bousquet | Introduction to statistical learning theory[END_REF] to provide bounds on the error made when working with samples instead of original trajectories.

Organization of the paper.. In Section 2 we provide the context and background on optimal control and occupation measures. In Section 3, we present our characterization of solutions to the inverse optimal control problem and illustrate how it allows to further discuss about the set of solutions and links with the direct optimal control problem. Numerical approximations via polynomials and statistical approximations via finite samples are provided and discussed in Section 4. The resulting numerical scheme (with proven strong theoretical guarantees) can be implemented with off-theshelf software on a standard computer. Finally, Section 5 describes numerical results on academic examples.

Preliminaries.

2.1. Notations. If A is a compact subset of a finite-dimensional Euclidean space, let C(A) resp. C 1 (A) denote the set of continuous resp. continuously differentiable functions from A to R. Let M(A) denote the space of Borel measures on A, the topological dual of C(A) with duality bracket denoted by ., . , i.e. µ, f = A f (x)dµ(x) is the integration on A of a function f ∈ C(A) with respect to a measure µ ∈ M(A). Let M + (A) resp. C + (A) denote the cone of non-negative Borel measures resp. non-negative continuous functions on A. The support of a measure µ ∈ M + (A) is denoted by spt µ. An element µ ∈ M + (A) such that µ, 1 = 1 is called a probability measure. Let δ x denote the Dirac measure concentrated on x and let I(e) denote the indicator function of an event e, equal to 1 if e is true, and 0 otherwise.

Let X ⊆ R d X denote the state space and U ⊆ R d U denote the control space which are supposed to be compact subsets of Euclidean spaces. System dynamics are given by a continuously differentiable vector field f ∈ C 1 (X × U ) d X . Terminal state constraints are modeled by a set X T ⊂ X which is also given. Let B n denote the unit ball of the Euclidean norm in R n , and let ∂S denote the boundary of set S in the Euclidean space. Let R[z] denote the set of multivariate polynomials with real coefficients with variables z and let R k [z] denote the set of such polynomials with degree at most k. For a polynomial p ∈ R k [z], we denote by p 1 the sum of the absolute values of the coefficients of p when expanded in the monomial basis.

2.2. Context: free terminal time optimal control. We consider direct optimal control problems of the form:

v 0 (z) := inf u,T T 0 l 0 (x(t), u(t))dt s.t. ẋ(t) = f (x(t), u(t)), x(t) ∈ X, u(t) ∈ U, t ∈ [0, T ], x(0) = z, x(T ) ∈ X T , T ∈ [0, T M ] (ocp 0 )
with Lagrangian l 0 ∈ C(X × U ) and free final time T with a given upper bound T M which ensures that the value function v 0 is bounded below. Dynamics f are given, as well as the sets X, U and X T ⊂ X. We assume that a set X 0 ⊂ X is given such that the following assumption is satisfied:

Assumption 1. For all initial conditions z ∈ X 0 , problem (ocp 0 ) is feasible.

Occupation measures.

In this section we describe how to construct an occupation measure from a feasible trajectory of (ocp 0 ) and then from a set of such trajectories. The content of this section was already described in the litterature (see for example [START_REF] Lasserre | Nonlinear optimal control via occupation measures and LMI relaxations[END_REF][START_REF] Henrion | Convex computation of the region of attraction of polynomial control systems[END_REF][START_REF] Henrion | Optimization on linear matrix inequalities for polynomial systems control[END_REF]) and we include these notions here for completeness. Let z ∈ X 0 be an initial point. We use Assumption 1 to fix a trajectory starting from z. That is, a terminal time T z ∈ R + , a measurable control u z : [0, T z ] → U and an absolutely continuous trajectory

x z : [0, T z ] → X such that ẋz (t) = f (x z (t), u z (t)), x z (0) = z, x z (T z ) ∈ X T .
(2.1)

The occupation measure of the corresponding trajectory is denoted by µ z and is defined by

µ z (A × B) := Tz 0 I(x z (t) ∈ A, u z (t) ∈ B)dt (2.2)
for every Borel sets A ⊂ X and B ⊂ U . We now turn to the construction of occupation measure and terminal measure of a set of trajectories by taking a measurable combination of occupation measures of single trajectories. Consider a probability measure µ 0 ∈ M + (X 0 ) and an upper bound on terminal time T M . Thanks to Assumption 1, for each z ∈ spt µ 0 , we fix a terminal time T z ∈ [0, T M ], a measurable control u z : [0, T z ] → U and an absolutely continuous trajectory x z : [0, T z ] → X such that (2.1) holds. That is, for each z ∈ spt µ 0 , we have an occupation measure µ z as described in (2.2). The occupation measure µ ∈ M + (X × U ) and terminal measure µ T ∈ M + (X T ) of the set of trajectories {x z (t)} z∈spt µ0,t∈[0,Tz] are then defined by tacking a convex combination of each µ z according to µ 0 (see also [START_REF] Henrion | Optimization on linear matrix inequalities for polynomial systems control[END_REF]Chapter 5] and [START_REF] Henrion | Convex computation of the region of attraction of polynomial control systems[END_REF]Section 3]). We obtain the following definition:

µ(A × B) := X0 µ z (A, B)µ 0 (dz), = X0 Tz 0 I(x z (t) ∈ A, u z (t) ∈ B)dt µ 0 (dz), µ T (A) := X0 I(x z (T z ) ∈ A) µ 0 (dz), (2.3)
for every Borel sets A ⊂ X and B ⊂ U . With the previous definition,

µ, l = X0 Tz 0 l(x z (t), u z (t))dt µ 0 (dz), ∀l ∈ C(X × U ).
In particular

µ(X × U ) = µ, 1 = X0 T z µ 0 (dz). Furthermore for every v ∈ C 1 (X), µ, grad v • f = X0 Tz 0 grad v(x z (t)) • f (x z (t), u z (t)) dt µ 0 (dz) = X0 (v(x z (T z )) -v(x z (0))) µ 0 (dz) = µ T , v -µ 0 , v , (2.4)
where "grad" denotes the gradient vector of first order derivatives of v, and the "dot" denotes the inner product between vectors. Equation (2.4) is known as Liouville's equation and is also written as

divf µ + µ T = µ 0 , (2.5)
where the divergence is to be interpreted in the weak sense and a change of sign comes from integration by part. As we have seen, occupation and terminal measures as defined in (2.3) satisfy the Liouville equation (2.5). This motivate the following broader definition. Definition 2.1. A general occupation measure is a measure that satisfies Liouville's equation (2.5), for some terminal measure µ T ∈ M + (X T ), in the weak sense described in (2.4).

We have seen in this section how to construct an occupation measure from a set of feasible trajectories of (ocp). However the set of all occupation measures is in general much bigger than the set of measures arising in this way.

2.4. Input of the inverse optimal control problem. For inverse optimal control, we suppose that the trajectories are given. Moreover, the Liouville equation and positivity constraints are sufficient to develop all the aspects of our analysis of inverse optimality.

Therefore, independently of how it is constructed, the input data of our inverse control problem is a general occupation measure as given by Definition 2.1.

This restriction is made without loss of generality regarding classical trajectories because, from the construction in (2.3), we consider an input set that contains all of them. All the results will in particular apply to situations when the occupation measure is a superposition of classical trajectories as described in (2.3). The results will also hold if this is not the case and the input measure involves generalized control. Finally and most importantly, this construction allows to formally treat cases for which we are given a possibly uncountably infinite number of trajectories as input data and is therefore much more general than considering one or a few classical trajectories.

2.5. Direct optimal control. Using the formalism of occupation measures, given a continuous Lagrangian l, an initial measure µ 0 and a maximal terminal time T M , we consider direct optimal control problems of the form

p * (µ 0 ) := inf µ,µ T µ, l s.t. divf µ + µ T = µ 0 , µ, 1 ≤ T M , µ ∈ M + (X × U ), µ T ∈ M + (X T ). (ocp) Definition 2.2 (OCP). OCP(l, µ 0 , T M ) is the set of measures (µ, µ T ) solving problem (ocp).
Note that by Lemma 2.3 and Assumption 1, set OCP(l, µ 0 , T ) is not empty. The link between problems (ocp 0 ) and (ocp) is far from trivial. . It is possible to construct problems for which measures considered in problem (ocp) do not arise in this way which may introduce spurious minimizers which are far from classical trajectories of problem (ocp 0 ), see for example [START_REF] Henrion | Convex computation of the region of attraction of polynomial control systems[END_REF]Appendix C]. These problems are usually overly constrained and not physically relevant, and in most practical settings, we have

p * (δ z ) = v 0 (z) ∀ z ∈ spt µ 0 ,
which we could see as an assumption on the inverse problem data. In this constrained setting, sufficient conditions for this property to hold are those that ensure the applicability of the Filippov-Ważewski Theorem, see [START_REF] Frankowska | Filippov's and Filippov-Wazewski's theorems on closed domains[END_REF] and the discussion around [17, Assumption I], [19, Assumption 2], [20, Assumption 1]. Under such sufficient conditions, it can be shown using [START_REF] Vinter | Convex duality and nonlinear optimal control[END_REF]Theorem 2.3] that the equality holds. However, as we argued in the introduction, the link between (ocp 0 ) and (ocp) is much less problematic when considering inverse optimality. The main reason is that we consider that the input of the inverse problem is a measure, which is therefore given and fixed. It could arise as in (2.3) but not necessarily (see Figure 1.1). We would like to emphasize the following:

• if the input occupation measure does not satisfy (2.3), then it does not make sense to consider (ocp 0 ) as a basis for inverse optimality since the input of the problem itself is more general than the classical controls considered in (ocp 0 ). In this case, it is more relevant to focus on (ocp) only. • if the input occupation satisfy (2.3), then, the analysis is still valid. In this case, since the input of the inverse problem involves classical controls, the question of the link between (ocp 0 ) and (ocp) is a real issue for direct optimal control. But in the context of inverse optimality, a partial answer is given a posteriori by Corollary 3.3. It is shown that, even in this case, considering (ocp) as a basis for inverse optimality does not allow to identify Lagrangians for which there is a gap between (ocp 0 ) and (ocp) for all considered initial conditions in spt µ 0 , except for a µ 0 -negligible subset.

For these reasons we adopt the following convention

All our analysis refers to direct control problems of the form of (ocp).

and the link with (ocp 0 ) (when it makes sense) will be a posteriori justified by Corollary 3.3: The corresponding conic dual can be written as

d * (µ 0 ) := sup v,w µ 0 , v -wT M s.t. l + w + grad v • f ∈ C + (X × U ) -v ∈ C + (X T ), w ≥ 0, v ∈ C 1 (X), w ∈ R. (hjb)
The first two constraints l

+ w + grad v • f ∈ C + (X × U ) and -v ∈ C + (X T ) of (hjb)
are relaxations of the well-known Hamilton-Jacobi-Bellman (HJB) sufficient condition of optimality [START_REF] Athans | Optimal control. An introduction to the theory and its applications[END_REF][START_REF] Bardi | Optimal control and viscosity solutions of Hamilton-Jacobi-Bellman equations[END_REF]. Conic duality provides the following link between the problems (ocp) and (hjb). Lemma 2.3. The infimum in (ocp) is attained and there exists a maximizing sequence in (hjb). In addition, for any feasible primal pair (µ, µ T ) and any sequence of dual variables v k ∈ C 1 (X) and w k ∈ R, k ∈ N, the following assertions are equivalent

• (µ, µ T ) is optimal for (ocp) and (v k , w k ) k∈N is a maximizing sequence for (hjb); • strong duality:

µ 0 , v k -w k T M -→ k→∞ µ, l ; (2.6)
• complementarity:

w k ( µ, 1 -T M ) -→ k→∞ 0, µ, l + w k + grad v k • f -→ k→∞ 0, µ T , v k -→ k→∞ 0.
(2.7)

Proof. We only sketch the proof here, for more details see [START_REF] Lasserre | Nonlinear optimal control via occupation measures and LMI relaxations[END_REF]. Observe that (ocp) is feasible thanks to Assumption 1 and (hjb) is feasible with w = max(-min X×U l, 0) and v = 0. Moreover, the cone

{(divf µ + µ T , µ, 1 , µ, l 0 ) : µ ∈ M + (X × U ), µ T ∈ M + (X T )} is closed for the weak topology σ(M + (X) × R 2 , C + (X) × R 2 ) (
by using Banach-Alaoglu's Theorem). Therefore there is no duality gap between (ocp) and (hjb) and the optimum is attained in the primal, see e.g. [7, Theorem IV.7.2]. Condition (2.6) is just a reformulation of strong duality in this context. Equivalence with (2.7) follows by noticing that for any primal feasible pair (µ, µ T ) and dual feasible pair (v, w),

µ, l = µ, l -divf µ + µ T -µ 0 , v ≥ µ, l -divf µ + µ T -µ 0 , v + w( µ, 1 -T M ) = µ, l + w + grad v • f -µ T , v + µ 0 , v -wT M ≥ µ 0 , v -wT M .
Remark 1. If the Lagrangian l is strictly positive on X × U , then Lemma 2.3 holds without the constraint µ, 1 ≤ T M and without the dual variable w.

3. Inverse optimal control. Given a "set" of trajectories and model constraints, the inverse problem of optimal control consists of finding a Lagrangian for which the trajectories are optimal. Thanks to the framework exposed in the previous section, it is now easy to define what is a solution to the inverse optimal control problem.

• Firstly, the "set" of trajectories will be represented by measures satisfying Liouville equation (2.5) which are part of the data of the inverse problem.

• Secondly, a Lagrangian l solution to the inverse problem is a continuous function such that (µ, µ T ) ∈ OCP(l, µ 0 , T ) for some T such that OCP(l, µ 0 , T ) is feasible.

In this section, we propose a rigorous definition of inverse optimality and prove an equivalence result between direct and inverse optimality. To do so, we use Lemma 2.3 which ensures that OCP(l, µ 0 , T ) is non empty as long as T ≤ T M . Furthermore, it provides a certificate of (sub)optimality.

3.1.

What is a solution to the inverse optimal control problem?. We can now formally define what is meant by a solution to the inverse optimal control problem:

Definition 3.1 (IOCP and IOCP ). For > 0, given measures µ ∈ M + (C × U ) and µ T ∈ M + (X T ) such that divf µ + µ T ∈ M + (X 0 ), denote by IOCP (µ, µ T ) the set of -optimal solutions to the inverse optimal control problem, namely the set of functions l ∈ C(X × U ) such that there exists a function v ∈ C 1 (X) satisfying

µ, l + grad v • f ≤ , l + grad v • f + ∈ C + (X × U ), µ T , v ≥ -, -v ∈ C + (X T ).
Then the set IOCP(µ, µ T ) of solutions to the inverse optimal control is defined by:

IOCP(µ, µ T ) := {l ∈ C(X × U ) : l ∈ IOCP (µ, µ T ) ∀ > 0 }.
Intuitively, Definition 3.1 states that we can find differentiable suboptimality certificate for any arbitrary precision (see in Remark 2). In addition, the positivity constraint on l + grad v • f + ensures that these certificates provide lower bounds on the value of the direct problem (ocp 0 ) for arbitrary initial conditions, even not in spt µ 0 . The main motivation behind this definition of inverse optimality is the following:

Theorem 3.2. Given µ ∈ M + (C × U ) and µ T ∈ M + (X T ), the set IOCP(µ, µ T
) is a convex cone, closed for the supremum norm. Moreover, the following two assertions are equivalent:

• l ∈ IOCP(µ, µ T ), divf µ + µ T = µ 0 ∈ M + (X 0 ); • ∃ T > µ, 1 , (µ, µ T ) ∈ OCP(l, µ 0 , T ). Proof.
Convexity follows from convexity of the constraints of Definition 3.1. There exists a constant K such that for any pair (l 0 , v 0 ) that satisfies constraints of Definition 3.1 for a certain > 0, then it holds for any Lagrangian l that µ, l

+ grad v 0 • f ≤ + K l -l 0 ∞ and l + grad v 0 • f ≥ --K l -l 0 ∞ on X × U , which is sufficient to prove closedness.
For the first implication, suppose that l ∈ IOCP(µ, µ T ) and divf µ + µ T = µ 0 . Then for any T > µ, 1 , the pair (µ, µ T ) is feasible for OCP(l, µ 0 , T ). Lemma 2.3 holds and the definition of IOCP(µ, µ T ) allows to construct a dual sequence that is feasible for OCP(l, µ 0 , T ) and that satisfies the complementarity condition with the pair (µ, µ T ).

We now turn to the last implication. Suppose that (µ, µ T ) ∈ OCP(l, µ 0 , T ) and µ, 1 < T . In particular, (µ, µ T ) is feasible for OCP(l, µ 0 , T ) and Lemma 2.3 holds. Consider the dual maximizing sequence (v k , w k ) k∈N given by Lemma 2.3. Complementarity ensures that lim k→∞ w k = 0. Furthermore, it holds that lim k→∞ µ, l + grad

v k • f = 0, l + grad v k • f ≥ -w k on X × U , lim k→∞ µ T , v k = 0 and v k ≤ 0 on X T which shows that l ∈ IOCP(µ, µ T ).
Remark 2. Another motivation behind Definition 3.1 of inverse optimality is the following. Suppose that l ∈ IOCP (µ, µ T ), then for any T ≥ µ, 1 , (µ, µ T ) is close to optimal for the problem OCP(l, µ 0 , T ). Indeed, suppose that (μ, μT ) ∈ OCP(l, µ 0 , T ). Then there exists v such that

µ, l + grad v • f = µ, l + µ T -µ 0 , v ≤ μ, l + grad v • f = μ, l + μT -µ 0 , v ≥ -T ,
and -µ T , v ≤ as well as μT , v ≤ . In addition, µ, l ≤ µ 0 , v + 2 and μ, l ≥ µ 0 , v -T . Therefore µ, l ≥ μ, l ≥ µ, l -(T + 2) .

Remark 3. At first sight the introduction of T in Theorem 3.2 may look artificial whereas in fact it carries important information. The second part in the equivalence states that (µ, µ T ) is a solution to some direct problem and does not saturate one of the constraints. This allows to avoid direct problems for which, for any value of T , any solution would saturate the constraint on the mass of the occupation measure; for example this happens in direct problems with free terminal time tending to infinity. Such problems should be avoided since then an occupation measure with finite mass cannot be optimal. Given a Lagrangian l, there is no guarantee that there exists a triplet (µ, µ T , T ) which satisfies the second point of Theorem 3.2. However, checking that a Lagrangian l meets our criterion for inverse optimality ensures that this is the case.

An interesting corollary is that if the input of the optimal control is given by classical trajectories, then inverse optimality ensures that the value of (ocp 0 ) is attained by classical trajectories for almost all the initial values considered. This leaves aside many of the technical issues when working with classical trajectories for direct optimal control. Corollary 3.3. If l ∈ IOCP(µ, µ T ) and (µ, µ 0 , µ T ) is a superposition of classical trajectories as defined in equation 2.3 in Section 2.3, then µ 0 -a.a. (almost all) of these trajectories must be optimal for the corresponding direct problem. In particular, v 0 (z) given by (ocp 0 ) is attained and there is no relaxation gap between (ocp 0 ) with initial condition z and (ocp) with initial measure δ z for µ 0 -a.a. initial conditions z in sptµ 0 .

As a consequence, the focus on (ocp) instead of (ocp 0 ) in Section 2.5 is a posteriori justified by Corollary 3.3. The question of absence of such a gap for initial conditions z ∈ spt µ 0 cannot be treated by this approach. This question is much less relevant for inverse optimality since it does not involve initial conditions that are related to input data of the inverse problem.

3.2. Applications to inverse optimality. We claim that Definition 3.1 is a powerful tool to analyze inverse optimality in the context of optimal control. To go beyond Theorem 3.2, we next describe results and comments that stem from Definition 3.1 of inverse optimality.

3.2.1.

How big is the space of solutions to the inverse problem?. Theorem 3.2 justifies the idea that if trajectories realize the minimum of some optimal control process then the corresponding Lagrangian meets our criterion. This requirement is necessary for any "inverse problem" (and not only for inverse optimal control). However in general there could be many candidate solutions as illustrated in this section. In what follows, we assume that the triplet (µ, µ 0 , µ T ) satisfies Liouville's equation (2.5).

Conserved values.. Suppose that there exists a function g ∈ C(X × U ) such that g(x, u) = 0 for all (u, x) ∈ spt µ. Then g 2 ∈ IOCP(µ, µ T ). In practical examples there might be many such conserved values. For instance this is the case when x or u or both remain on a manifold or when there exists a continuous mapping x → u(x).

Total variations.. Consider any function g ∈ C 1 (X). All Lagrangians of the form l = grad g • f belong to IOCP(µ, µ T ), independently of (µ, µ 0 , µ T ).

Convex conic combinations and uniform limits of solutions.. As stated in Theorem 3.2, the set of solutions to the inverse problem is a convex cone, closed for the supremum norm. For example, let g ∈ C 1 (X) and consider a Lagrangian l ∈ C(X ×U ). Then OCP(l, µ 0 , T ) = OCP(l + grad g • f, µ 0 , T ) for every T > 0 and therefore both Lagrangians l and l + grad g • f are solutions to the inverse optimal control problem.

All the above examples illustrate that many solutions to the inverse problem may exist. Although these solutions are valid from a theoretical point of view, they do not correspond to what is commonly expected from a solution. Indeed, they do not arise from an optimal physical process that would have generated trajectories, but rather from mathematical artifacts.

3.2.2. How does the direct problem affect the space of solutions to the inverse problem?. Intuitively, the more information is contained in (µ, µ T , µ 0 ), the smaller is the space of solutions to the inverse problem. We next discuss two factors that impact the size of IOCP(µ, µ T ).

Direct problem constraints.. Denote by K 1 (resp. K 2 ) the feasible set of problem (ocp) and by IOCP 1 (µ, µ T ) (resp. IOCP 2 (µ, µ T )) the set of solutions to the inverse problem (as described in Definition 3.1) when the state, control and dynamical constraints are given by (X 1 , U 1 , f 1 ) (resp. (X 2 , U 2 , f 2 )).

If K 1 ⊂ K 2 then IOCP 2 (µ, µ T ) ⊂ IOCP 1 (µ, µ T ). In other words, there is a kind of duality between the space of feasible solutions for the direct problem and the space of solutions to the inverse problem. An extreme instance is when the feasible space of the direct problem is a singleton (f does not depend on the control u), in which case any Lagrangian is a solution to the inverse problem.

Range of the occupation measure.. Suppose that (µ, µ 0 , µ

T ) = (µ 1 , µ 1 0 , µ 1 T ) + (µ 2 , µ 2 0 , µ 2 T ) where divf µ 1 + µ 1 T = µ 1 0 and divf µ 2 + µ 2 T = µ 2 0 . Then IOCP(µ, µ T ) ⊂ IOCP(µ 1 , µ 1
T ). As a consequence, maximizing the support of the initial measure µ 0 reduces the space of solutions to the inverse problem. When the occupation measure µ is a superposition of trajectories as detailed in Section 2.3, the larger is the "space" occupied by trajectories, the smaller is the space of potential solutions to the inverse problem.

A toy example of quantitative well-definedness analysis.

To illustrate the proposed framework we consider a simple uni-dimensional example. We emphasize that his example is very simple in the sense that the direct problem is easy. However, inspecting the solution of the inverse problem leads to non trivial behaviors. Let X = [-1, 1] with X T = {0} and let f (x, u) = u with U = [-1, 1]. Consider the family of Lagrangians F = {l α : u → 1 + α 2 u 2 , α ≥ 0}. Suppose that we are given a triplet (µ, µ 0 , µ T ) which consists of a superposition of trajectories as described in Section 2.3. We wish to find a candidate Lagrangian in the family F. Then we have the following alternatives.

1 F ∩ IOCP(µ, µ T ) = ∅. 2 F ∩ IOCP(µ, µ T ) = F. 3 F ∩ IOCP(µ, µ T ) is a singleton, {l α }, α > 2. 4 F ∩ IOCP(µ, µ T ) = {l α , 0 ≤ α ≤ 2} .
We should comment on case 1 latter. If the support of µ is empty, which means that µ 0 = µ T = δ x , then we are in case 2. Assume now that the support of µ is non-empty and we are not in case 1. Then there exists α ≥ 0 such that l α ∈ IOCP(µ, µ T ). Consider a sequence of decreasing positive numbers k → 0 and the corresponding certificates functions v k that allow to verify that l α ∈ IOCP k (µ, µ T ). Since µ T = δ 0 , we may assume (up to an addition) that v k (0) = 0 which simplifies the problem. In addition, one must have 1 + α 2 u 2 + v k (x)u → 0, µ almost every where (recall that the support of µ is non empty and this concerns a non empty subset of X and U ). Furthermore, for any x, one must have 1

+ α 2 u 2 + v k (x)u + k = α 2 u + v (x) √ 2α 2 + 1 + k -v (x) 2 2α ≥ 0, for u ∈ [-1, 1]. • Suppose that α > 2. Then for any k, |v k (x)| ≤ 1 + α 2 + k .
Since k goes to 0 and α > 2, for k sufficiently large, has and therefore |u| = 1, µ almost everywhere. In this case, it is easy to construct alternative sequences ṽk = 2+ α 2+ α v k for 0 ≤ α ≤ 2 to show that l α is also a member of IOCP(µ, µ T ). To conclude, we are in case 1 when the trajectories that generate µ are not optimal with respect to any Lagrangian in F, in particular when |u| is not µ almost everywhere constant. If this is not the case and µ is not degenerate, we have a unique solution or a set of solutions depending on µ and its relation with the constraint on u.

|v k (x)| α ≤ 1. Taking u = - v k (x) α gives 1 + α 2 u 2 + v k (x)u + k ≥ 1 + k -v (x) 2 2α ≥ 0. It must hold µ almost everywhere that 1 + k -v (x) 2 2α → 0 and |u| = 2 α . • Suppose that 0 ≤ α ≤ 2. It must hold µ almost every where that v k (x)u → -1 -α 2 u 2 .

Practical inverse control.

As discussed in Section 3.2.1, the space of solutions to the inverse problem can be very large. Many of these solutions are of little interest for practitioners because they lack some physical meaning. However, from a formal point of view "valid" solutions exist and ideally they should be the only solutions of a practical inverse optimal control problem to be defined.

One may invoke some heuristics to reduce the space of solutions and to enforce prior knowledge in the treatment of the inverse problem. This is commonly achieved by imposing constraints on the candidate Lagrangian solution. Such heuristics include :

• restricting the dependence on certain variables;

• shape conditions (e.g., convexity);

• conic constraints such as positivity;

• parametric constraints (e.g., considering a finite dimensional family of candi-date Lagrangians); • constraints relating the dependence between the candidate Lagrangian and the corresponding value function. Notice that Definition 3.1 refers to the large class of continuous Lagrangians with conic constraints. From a theoretical perspective, this allows to characterize inverse optimality in full generality. However this is not amenable to numerical computation yet and so we also describe tractable numerical approximations in the context of inverse optimality.

Finally, according to Definition 3.1, the input of the inverse problem is an occupation measure. Again, this is a convenient tool for theoretical purposes but in most practical cases such an occupation measure is not available. In fact, roughly speaking, only some realizations of an experiment are available and these realizations form a data set which is an approximation of an hypothetical occupation measure. Therefore in practice the input of the inverse problem is only an approximation of an ideal input, and correctness of this approximation is justified under certain experimental assumptions at the end of this section.

Normalization.

The trivial Lagrangian is solution to the inverse problem independently of the input occupation measures. As we have seen in Section 3.2.1, total variations share the same property. Even though these are solutions to the inverse problem, it is important to avoid them in practice because they do not depend on the input occupation measure and therefore carry no information about it. As illustrated in Example 3.2.3, one way to avoid these spurious solutions is to consider only very restricted families of Lagrangians that cannot contain such solutions. This might be quite restrictive in practice and therefore we provide an alternative. We need the following assumption:

Assumption 2 (Finite time controllability). There exists T > 0 and a compact set X ⊂ X with smooth boundary and µ 0 ( X) > 0, such that for any x 1 , x 2 ∈ X, there exists s ∈ [0, T ], a bounded function u : [0, s] → U and an absolutely continuous trajectory x : [0, s] → X such that x(0) = x 1 , x(s) = x 2 and ẋ(t) = f (x(t), u(t)) for all t ∈ [0, s].

Under assumption 2 we have the following result. Proposition 4.1. If in Definition 3.1 one includes the normalization

1 - X×U l + grad v • f = then 0 ∈ IOCP(µ, µ T ).
Proof. This is due to the following contradiction. Suppose that 0 ∈ IOCP(µ, µ T ). Choose a decreasing sequence k → 0 as k → ∞, and construct a sequence v k of differentiable functions that satisfy conditions of Definition 3.1 for the chosen k . Then

µ 0 , v k = µ T , v k -µ, grad v k • f = O( k ). Because of Assumption 1, v k (x 0 ) ≤ T M k for every x 0 ∈ X 0 . Therefore, by integration, max( X v k µ 0 , X\ X v k µ 0 ) ≤ O( k ). Furthermore, X v k µ 0 + X\ X v k µ 0 = µ 0 , v k = O( k ). Finally, max X v k µ 0 , - X v k µ 0 ≤ O( k )
and lim k→∞ X v k = 0. In addition, by Assumption 2, we also have

|v k (x 1 ) -v k (x 2 )| ≤ T k for every x 1 , x 2 ∈ X. Since µ 0 ( X) > 0, this implies that v k → 0 uniformly on X, as k → ∞.
In addition, X×U grad v k • f → 1 as k → ∞. Next, with the polynomial vector field

x → g(x) = U f (x, u)du, Stokes' Theorem yields X div(v k (x)g(x)) = X div(g(x))v k (x) + grad v k (x) • g(x) = X div(g(x))v k (x) + X×U grad v k (x) • f (x, u) = ∂ X v k (x)g(x) • n(x)
where n(x) is the outward pointing normal to the boundary at x. Because of the uniform convergence of v k → 0 on X as k → ∞, and boundedness of g and div g on X, the left-hand side converges to 1 while the right-hand side converges to 0, which is a contradiction. Remark 4. The conditions on X may be relaxed. Indeed, the only important point is to be able to apply Stokes's Theorem. In particular, the set X could be a box or an open set whose boundary does not have too many non-smooth points, see e.g. [START_REF] Whitney | Geometric integration theory[END_REF]Theorem III.14A].

Remark 5. The normalization given in Proposition 4.1 obviously ensures that Lagrangians in the form of a total variation are excluded. Assumption 2 may look very strong regarding the result of Proposition 4.1. However the next example shows that it cannot be excluded.

Example 1. Consider the direct control problem with X = U = [0, 1], X T = {1}, T M = 1 and f (x, u) = u. These data are obviously not compatible with Assumption 2. Choose µ 0 (dx) = dx and l 0 (x, u) = 1 so that the couple µ(dx, du) := xdxδ 1 (du) and µ T (dx) := δ 1 (dx), solves the problem

p * (µ 0 ) := inf µ,µ T µ, 1 s.t. divf µ + µ T = µ 0 , µ, 1 ≤ T M , µ ∈ M + (X × U ), µ T ∈ M + (X T ). (4.1) Indeed, for any differentiable v, divf µ, v = 1 T -v f µ(dx) = 1 T vdx -v(1) = µ 0 , v -µ T , v . Furthermore, the function x → v * (x) := 1 -x ensures that (µ, µ T ) is an optimal solution. Indeed µ, 1 = µ 0 , v * . Consider a sequence of differentiable functions (v k ), k ∈ N, such that v k (x) = 0 for x ≥ 1 k and v k (x) = -(kx -1) 2 other- wise. For X = [0, 1], we have 0 ≤ µ, v k f = 1 k T kx -k 2 x 2 dx = 1 6k → 0, v k f ≥ 0, µ T , v k = 0, v k (1) = 0, X×U v k f = v k (1) -v k (0) = 1.
Therefore even if we enforce the normalization X×U l + v f = 1, we cannot prevent the trivial Lagrangian l = 0 from being an optimal solution to the inverse problem. Remark 6. Regarding Proposition 4.1, one could argue that simple linear constraints such as l(0) = 1 or conic constraints such as l > 0 would be sufficient to avoid the trivial Lagrangian. However, this does not allow to avoid total variations which are equivalent to the trivial Lagrangian in terms of solutions to the direct problem.

Denote by NIOCP the subset of IOCP with the normalization constraint of Proposition 4.1 added to the constraints of Definition 3.1. One important feature of this normalization is that it can be thought of as a way to intersect the cone of solutions to the inverse problem with an affine subspace. Therefore NIOCP is still closed and convex. Furthermore, if we restrict the set of candidate Lagrangians to be finitedimensional, one may look for minimum norm-like solutions which will prove to be useful in numerical experiments. Indeed, NIOCP being closed and convex and all norms being equivalent in finite dimensions, optimization problems over NIOCP, if bounded, have an optimal solution. Finally, as 0 ∈ NIOCP, one may minimize any norm-like function to enforce specific prior structure and avoid the trivial Lagrangian.

Polynomial approximation.

Until now, all the results that we have presented involve continuous and differentiable functions, in full generality. However, for practical computation one has to approximate such functions and of course, polynomials are obvious natural candidates. But in our context they are also of particular interest for mainly three reasons:

• for fixed degree, polynomials belong to finite-dimensional spaces and are therefore amenable to computation; • when varying the degree, the class of polynomials is rich enough to approximate a wide class of functions; • Positivity Certificates from real algebraic geometry allow to express positivity constraints in a computationally tractable way. From now on, we make the following assumption: Assumption 3. f is a polynomial and X, X T and U are basic semi-algebraic sets. As proposed in Definition 3.1, checking that a polynomial is a solution of the inverse problem involves the construction of a sequence of continuously differentiable functions. These functions can also be approximated by polynomials. In this section we describe some tools required for such an approximation and we also prove the correctness of the approximations.

Let g 1 , . . . , g m ∈ R[z] be polynomials in the variable z ∈ R n and consider the basic semi-algebraic set

G = {z ∈ R n : g i (z) ≥ 0, i = 1, . . . , m}.
Let g 0 = 1 and let Σ 2 ⊂ R[z] denotes the set of sums-of-squares (SOS) polynomials, i.e., p ∈ Σ 2 if it can be written as a sum of squares of other polynomials. Definition 4.2. Let Q k (G) denote the convex cone of polynomials that can be written as

p = m i=0 s i g i , s i ∈ Σ 2 , i = 0, 1, . . . , m,
where the degree of s i g i , i = 0, 1, . . . m, is at most 2k. If p ∈ Q k (G) we say that p has a Putinar positivity certificate. It is immediate to check that any element of Q k (G) is non-negative on G. A remarkable property of such certificates is that a partial converse is true. Proposition 4.3 ([37]). Suppose that the polynomial super-level set {x :

g i (x) ≥ 0} is compact for some i = 1, . . . , m. If p > 0 on G then there exists k ≥ 0 such that p ∈ Q k (G).
Furthermore and importantly from a computational viewpoint, checking whether p ∈ Σ 2 reduces to checking whether a set of LMIs [START_REF] Lasserre | Moments, positive polynomials and their applications[END_REF] involving the coefficients of p In this section we first describe how the framework of occupation measures can formally describe the process of sampling trajectories and we justify the replacement of measures (µ, µ T ) by their empirical counterparts when considering empirical samples as input data for inverse control problems. In the context of polynomial certificates in Definition 4.4, this amounts to replacing the moments of the measures (µ, µ T ) by their empirical counterparts. Consider the probability measure µ 0 on X 0 and the measures (µ, µ T ) in equations (2.1) and (2.3). One way to interpret these measures is to consider the following random process:

• choose z ∈ spt µ 0 randomly, • choose t randomly uniformly on [0, T z ], • output ξ = (x z (t), u z (t)).
This defines a generative process for the random variable ξ. If the probability for an initial condition z to belong to a Borel set A ⊂ X 0 is given by

P[z ∈ A] =p 0 (A) := A T z µ 0 (dz) X0 T z µ 0 (dz)
, then the probability for a trajectory ξ to belong to a Borel hyperrectangle (A, B) ⊂ X × U is given by

P[ξ ∈ A × B] =p µ (A × B) := X0 Tz 0 I(x z (t) ∈ A)I(u z (t) ∈ B)dt µ 0 (dz) = µ(A × B) µ(X × U ) .
A statistical model for points (x i , u i ), i = 1, . . . , n, that are samples of trajectories, is to assume that we repeat the previous process n times, independently. In this case, we say that the database D = {(x i , u i )} n i=1 is made of independent realizations of a random variable with underlying distribution p µ . The process which generates the database being random, we write

{ξ i } n i=1 i.i.d. ∼ p µ .
to stress that all ξ i are independent and identically distributed (i.i.d.) according to p µ (they are independent copies of the same random variable). Similarly, µ T can be seen as the probability distribution describing the following random process:

• choose z ∈ spt µ 0 randomly according to p 0 , • output x z (T z ).
We now define what is an approximate solution to the inverse problem when the only information available about µ is a realization, D = {(x i , u i )} n i=1 , of a random process,

{ξ i } n i=1 i.i.d. ∼ p µ . Definition 4.6 (Sampled-IOCP ,k ). For > 0 and D = {(x i , u i )} n i=1 , let
Sampled -IOCP ,k be the set of polynomials l ∈ R[x, u] 2k such that :

1 n n i=1 l(x i , u i ) + grad v(x i ) • f (x i , u i ) ≤ , l + grad v • f ∈ Q k (X × U ), -v ∈ Q k (X T ), v + ∈ Q k (X T ), X×U l + grad v • f = 1.
for some polynomial v ∈ R[x] 2k . In other words, Sampled -IOCP ,k is the set of polynomial -optimal solutions (with degree bound k) of the sampled inverse optimal control problem., One has replaced µ by its empirical counterpart, added the normalization of Proposition 4.1 and simplified other conditions; see also Remark 7.

Importantly, membership in Sampled -IOCP ,k (D) can be tested by semi-definite programming. Using arguments from empirical processes and learning theory, one can quantify the price to pay for this discretization. Of course since we assume that the process that generates the data include some randomness, such a quantification holds probabilitically.

Proposition 4.7. Suppose that D = {(x i , u i )} n i=1 is a realization of the random process {ξ i } n i=1 iid ∼ p µ . Then there exist constants K 1 (X, U, k), K 2 (X, U, k) that only depend on (X, U ), and k ∈ N, such that for any l ∈ Sampled -IOCP ,k (D) and any 0 < δ < 1: l ∈ polyIOCP (n),k (µ, µ T ) with probability 1 -δ, where

(n) = + 1 √ n K 1 (X, U, k) + K 2 (X, U, k) ln 2 δ ,
and where the randomness comes from the realization of D. Proof. Apply Lemma A.1 of Appendix A to the polynomial l + grad v • f to get a bound on p µ , l + grad v • f . Remark 7. The conditions detailed in Definition 4.6 ensure that v ≤ 0 on X T and µ T , v ≥ -, for any terminal measure. This is done in order to avoid to deal with the terminal measure µ T but other alternatives are possible. For example, when X T is a single point or a simple algebraic set one may enforce (as we do in Section 5) v = 0 on X T instead. Another possibility is to replace µ T by its empirical counterpart. However in this case we need to provide a lower bound or add constraints on v to obtain finite sample bounds as described in Proposition 4.7. Proposition 4.7 mixes arguments from measure theory and conic optimization with arguments from empirical process and statistical learning theory. The implication of this result is that for a fixed degree, provided that the sample size is big enough, with high probability we do not loose much by approximating µ by an empirical sample.

• We would like to emphasize here that it is necessary to restrict the complexity of the class of functions in which the candidate Lagrangian is searched. Indeed otherwise for instance, for any fixed sample, the polynomial (x, u) → n i=1 ||x -x i || 2 ||u -u i || 2 belongs to Sampled -IOCP 0,4n (D), but this clearly does not give much insight on the original control problem! The degree of the polynomial candidates is one among many possible measures of complexity. Furthermore, although the constants are likely to be sub-optimal, they give a sense of how fast the degree of the polynomial approximation may grow with respect to the sample size in order to maintain accurate approximations of µ.

5. Numerical illustrations. Building on results of Section 4, we next provide illustrative numerical simulations. In order to fit in the framework of the previous section, we consider examples where f is a polynomial and X, U and X T are basic semi-algebraic sets. In addition, the input data of the inverse problem is given by a finite database: D = {(x i , u i )} i=1,...,n . In the sequel, the candidate Lagrangians satisfy Definition 4.6.

To compute such Lagrangians, the main idea is to solve an optimization problem with fixed k ∈ N, where:

• l ∈ R[x, u] 2k , v ∈ R[x]
2k and > 0 are the decision variables,

• is the criterion to minimize,

• Sampled -IOCP ,k (D) is the projection on (l, ) of the set of feasible solutions (l, , v). In addition, we also include a sparsity inducing term in the criterion that will prove to be useful in numerical experiments.

+ λ||l|| 1 s.t. 1 n n i=1 l(x i , u i ) + grad v(x i ) • f (x i , u i ) ≤ , l + grad v • f ∈ Q k (X × U ), v = 0 on X T , X×U l + grad v • f = 1. (iocp) where l ∈ R[x, u] 2k , v ∈ R[x] 2k
, is a real, λ > 0 (fixed) is a given regularization parameter, and ||.|| 1 denotes the 1 norm of a polynomial, i.e. the sum of absolute values of its coefficients when expanded in the monomial basis. The first constraints come from Definition 4.6 and the last affine constraint is meant to avoid the trivial solution; see Proposition 4.1. The 1 norm is not differentiable around sparse vectors (with entries equal to zero) and has the sparsity promoting role to bias solutions of the problem towards polynomial Lagrangian solutions with few nonzero coefficients. This regularization affects problem well-posedness and will prove to be essential in numerical experiments.

5.1.2. Numerical implementation. Linear constraints are easily expressed in term of polynomial coefficients. A classical lifting allows to express the 1 norm as a linear program: for x ∈ R n , ||x|| 1 = min s n i=1 s i subject to s i ≥ x i and s i ≥ -x i , for all i. The Putinar positivity certificates can be expressed as LMIs [START_REF] Lasserre | Moments, positive polynomials and their applications[END_REF] whose size depends on the degree bound k. We use the SOS module of the YALMIP toolbox [START_REF] Löfberg | Pre-and post-processing sum-of-squares programs in practice[END_REF] to manipulate and express polynomial constraints at a high level in MATLAB. Solution for the one dimensional minimum exit time problem, effect of the regularization parameter λ. See Section 5.2 for problem details and comments. The first column is the distribution of the error . It represents the value of l(x, u) + grad v • f (x, u) as a function of x where u = sign x is the optimal control. The second column is a representation of the value function v and the third column is a representation of its derivative for solutions of problem (iocp) with and without regularization. We take 100 points on the segment. Lagrangian l and value function v are both polynomials of degree 16.

small problems. As shown in the numerical results section, we could handle problems with 5 variables and degree 10 with a reasonable amount of time and memory. To handle larger size problems, specific heuristics and techniques beyond the scope of this paper must be implemented. 5.1.3. General setting. We consider several direct problems of the same form as (ocp 0 ). That is, we give ourselves compact basic semi-algebraic sets X, U , X T , the dynamics f , and a Lagrangian l 0 . We take known examples for which the (direct) optimal control law can be computed and try to vary their degree of difficulty. Given these optimal state-control trajectories, we generate randomly n data points D = {(x i , u i )} i=1...n according to the random process described in Section 4.6. For a given value of λ and k, we compute a solution l of problem (iocp). Then we measure how l is close to l 0 by computing the following quantity (in the monomial basis):

min α ||l 0 -αl|| 2 ||l 0 || 2 = 1 - l 0 , l 2 ||l 0 || 2 2 ||l|| 2 2 1 2
.

(5.1)

We also report the value of in program (iocp). A larger value of means less reliable numerical certificates; see Remark 2.

5.2. Illustration on a one-dimensional example. First consider the eikonal problem of minimum exit time from the unit ball in the one-dimensional case. The data of the problem are

X = U = B 1 , X T = ∂X, l 0 = 1, f = u.
The optimal law for this problem is u = signx and the value function is v 0 (x) = 1-|x|. We sample 100 points uniformly in X and solve problem (iocp). We compare the choices λ = 0 (no regularization) and λ = 1. Results are presented in Figure 5 which displays the distribution of the error , the estimated value function v as well as its first derivative. Despite the simplicity of the problem, it is quite representative of the difficulties that arise in the context of inverse optimality. The first difficulty is the size of the set of solutions to the inverse problem:

• given any symmetric differentiable concave function v vanishing on {-1, 1}, the pair (l = |, v) solves problem (iocp) with = 0; • any positive polynomial on X × U vanishing if |u| = 1 solves problem (iocp) with = 0; • any Lagrangian of the form v (x)u solves problem (iocp) with = 0; • any convex combination of solutions of the types mentioned above also solves problem (iocp).

Even though formally accurate, these solutions form a relatively large set and do not carry any physical meaning. In the absence of any additional form of prior knowledge, it is impossible to discriminate between these solutions and the one that we wish to recover, namely l 0 . This is illustrated in Figure 5.1 where the red line (λ = 0) displays an example of value function v obtained with a very low value of . This is a very good certificate that our database D is close to optimal for the corresponding Lagrangian l. However, the estimated Lagrangian is far from the original one, namely l 0 = 1. Moreover, the shape of the value function is quite uncommon. This motivates the use of prior knowledge to bias the solutions of problem (iocp) toward a certain set of solutions. We use the 1 norm which tends to promote Lagrangians with few non-zero coefficients.

When λ = 1, the sparsity inducing effect of 1 -norm regularization allows to recover the true Lagrangian (l 0 = 1) which is indeed sparse. The solution of problem (iocp) involves a polynomial function v which should in principle be close to the true value function v 0 (x) = 1 -|x|. The v function displayed in Figure 1 is close to v 0 . However, v 0 is not smooth around the origin and therefore its derivative is harder to approximate by polynomials around this point. Hence the value of the error is higher around the origin. that there is no gap between solutions of the direct and inverse problem (Theorem 3.2). To the best of our knowledge this is the first result of this kind. In addition, in principle the proposed methodology is applicable to practical problems where we only have access to sample trajectories. We have also proposed numerical and statistical approximation procedures from which solid theoretical guaranties can be obtained. Finally we have illustrated our results on relatively simple (but not trivial) numerical examples of modest size.

One of the most striking aspects of the inverse problem is its set of valid solutions. Indeed, even for the simplest problems it is difficult to discriminate between physically meaningful Lagrangians and spurious mathematical solutions. For this reason, formulating the inverse problem as a well-posed problem (in particular with a unique solution) requires the introduction of strong prior knowledge -sometimes arguably too restrictive -about the nature of the Lagrangian to be recovered; see for example [START_REF] Ajami | How humans fly[END_REF]. However the proposed formulation based on relaxed HJB-optimality conditions, allows to get intuitions about characteristics that affect well-posedness of the problem.

This work is to be seen as a first step toward a theoretical and practical framework for the resolution of inverse problems in a variety of contexts. Further aspects of the problem have to be investigated within this realm. First, we only deal with deterministic trajectories. For practical purposes it is essential to consider the effect of experimental noise, both from theoretical and practical perspectives, and to determine to which extent and how the problem can be solved in this more difficult context. Second, we have proposed a numerical scheme to approximate solutions and show that it is effective on academic examples of modest size. Experimental validation of such approximations should be carried out on real world examples of larger size. This involves a lot of data processing and fine tuning for each specific example. In this perspective, humanoid robotics provides an active and attractive field of application [START_REF] Arechavaleta | An optimality principle governing human walking[END_REF][START_REF] Mombaur | From human to humanoid locomotion-an inverse optimal control approach[END_REF].

Observe that K d (Z) is a subset of a finite dimensional space and that for p ∈ K d (Z), we have ||p|| Z 1 = Z |p| = 1. Since all norms are equivalent, K d (Z) is bounded in any given norm on polynomials, in particular, the supremum norm. Therefore, the quantity (A.1)

The left hand side depends on the random draw of the sample {z i } n i=1 , but the right hand side is deterministic. We use a standard symmetrization argument to bound the expectation in the right hand side. Using the definition of and n , the convexity of the supremum and Jensen's inequality, we have that where the notation z refers to any other sample S = {z i } n i=1 drawn from P and n is the corresponding empirical measure. The iid assumption allows to flip z i and z i in the expectation. Let ξ i be Rademacher variables, i.e. random variables which take values in {-1, 1}, each with probability one half. We have is also finite and attained. We have that sup

p∈K d (Z) 1 n c(p) • n i=1 ξ i v(z i ) ≤ M d c (Z) n n i=1 ξ i v(z i ) 2 = M d c (Z) n n i=1 n j=1 ξ i ξ j v(z i ) • v(z j ).
Moreover, E ξ n j=1 ξ i ξ j v(z i ) • v(z j ) = n i=1 v(z i ) 2 (E ξ [ξ i ξ j ] = I(i = j)). Therefore, using Jensen's inequality (with concavity of the square root), we obtain

E ξ sup p∈K 1 n c(p) • n i=1 ξ i v(x i ) ≤ M d c (X) n n i=1 ||v(x i )|| 2 2 ≤ M d c (X)M d v (X) √ n .
Putting things together, using inequalities (A.1), (A.2), (A.3), we have that with probability 1 -δ, it holds c, q ≤ c n , q + 2

M d c (X)M d v (X) √ n + M d ∞ (X) 1 2n ln 2 δ .

Figure 1 . 1 .

 11 Figure 1.1.Direct optimal and inverse optimal control flow chart. System description is given by the dynamics f , the state constraint set X and control constraint set U which are all assumed to be fixed. We emphasize that the Lagrangian and the occupation measure have symmetric roles for the direct and inverse problems. We also note that the output of the inverse problem is a Lagrangian. Solving the direct optimal control problem for new initial conditions is an important question but remains secondary regarding inverse optimality which is the focus of this work.

5. 1 .

 1 Numerical experiments. 5.1.1. Problem formulation. We consider the following optimization problem: inf l,v,

Figure 5 . 1 .

 51 Figure 5.1.Solution for the one dimensional minimum exit time problem, effect of the regularization parameter λ. See Section 5.2 for problem details and comments. The first column is the distribution of the error . It represents the value of l(x, u) + grad v • f (x, u) as a function of x where u = sign x is the optimal control. The second column is a representation of the value function v and the third column is a representation of its derivative for solutions of problem (iocp) with and without regularization. We take 100 points on the segment. Lagrangian l and value function v are both polynomials of degree 16.

  is finite. We have that for all i, and any z 1 , . . . , z n , z i and any p ∈ K d (Z)|F (z 1 , . . . , z i , . . . , z n ) -F (z 1 , . . . , z i , . . . , z n )| ≤ sup p∈K d (Z) 1 n |p(z i ) -p(z i )| ≤ M d ∞ (Z) n .Therefore McDiarmid's inequality of Lemma A.2 applies to functionF with α = M d ∞ (X) n, and, for any q ∈ K d (Z), with probability 1 -δ, it holds thatn , q ≤ sup p∈K d (Z) n , p ≤ E z sup p∈K d (Z) n , p + M d ∞

  E z sup p∈K d (Z) n , p = E z sup p∈K d (Z) E z [ n , p ]n , p (A.2) ≤ E z,z sup p∈K d (Z) nn , p = E z,z sup p∈K d (Z) i ) -p(z i )

ξ

  p(z i ) -p(z i )) ≤E z,z ,ξ sup p∈K d (Z) i p(z i ) =2E z,ξ sup p∈K d (Z) 1 n c(p) • n i=1 ξ i v(z i ) .The quantity on the right hand side is known as the Rademacher complexity of the function class K d (Z). Intuitively, it measures to which extent elements of a function class correlate with random noise in a worst case scenario. The function p → ||c(p)|| 2 is a norm on polynomials and since K d (Z) is bounded, the quantity M d c (Z) := sup p∈K d (Z) ||c(p)|| 2 is finite. Moreover, since Z is compact, the quantity M d v (Z) := sup z∈Z ||v(z)|| 2

  Error versus regularization parameter λ. Problems details are given in the main text. Estimation error (est) is given in (5.1). Epsilon error (eps) is the value of in program (iocp). Trajectory sample size: 20 for A, A' and B, 50 for C. Degrees of l and v are 4 and 10 respectively.

			A			A'			B			C
		1e+00
													Error
	error	1e-02											est eps
		1e-04
		1e-06	1e-03	1e+00	1e-06	1e-03	1e+00	1e-06	1e-03	1e+00	1e-06	1e-03	1e+00
						lambda
	Figure 5.2.
													.1

A semi-definite program is a finite-dimensional linear optimization problem over the cone of non-negative quadratic forms for which powerful primal-dual interior-point algorithms are available[START_REF] Vandenberghe | Semidefinite programming[END_REF].

Acknowledgments. This work was partly funded by an award of the Simone and Cino del Duca foundation of Institut de France, a grant of the Gaspard Monge program (PGMO) of the Fondation mathématique Jacques Hadamard. Most of this work was carried out during Edouard Pauwels' postdoctoral stay at LAAS-CNRS. The authors would like to thank Frédéric Jean, Jean-Paul Laumond, Nicolas Mansard and Ulysse Serres for fruitful discussions.

has a solution. Therefore a more precise definition of inverse optimality in the context of polynomial Lagrangians is as follows (compare with Definition 3.1).

Definition 4.4 (polyIOCP and polyIOCP ,k ). For > 0, given measures µ ∈ M + (C × U ) and µ T ∈ M + (X T ) such that div f µ + µ T = µ 0 ∈ M + (X 0 ), denote by polyIOCP ,k (µ, µ T ) the set of polynomials l ∈ R[x, u] 2k (i.e., of degree at most 2k) such that:

for some polynomial v ∈ R[x] 2k . Denote also by polyIOCP(µ, µ T ) the set of polynomial solutions to the inverse optimal control problem: That is, l ∈ polyIOCP(µ, µ T ) if for any > 0 there exists k( ) such that l ∈ polyIOCP ,k( ) .

In other words, polyIOCP ,k (µ, µ T ) is the set of polynomial -solutions with degree bound 2k, to the inverse optimal control problem. The advantage of the previous definition, is that, provided that one has the possibility to compute the linear functionals µ ∈ R[x, u] * and µ T ∈ R[x] * , checking whether l ∈ polyIOCP ,k (µ, µ T ) for k and given reduces to solving a convex LMI problem [START_REF] Lasserre | Moments, positive polynomials and their applications[END_REF]. Furthermore, under a compactness assumption, in the asymptotic regime this definition is equivalent to Definition 3.1.

Proposition 4.5 (Correctness of polynomial approximation). Suppose that one of the polynomials defining the basic semi-algebraic set X (resp. U ) has a compact super-level set. Then

Proof. The direct inclusion is trivial. For the reverse inclusion, suppose that l ∈ IOCP(µ, µ T )∩R[x, u]. Fix > 0, and take for v a certificate that l ∈ IOCP 3 (µ, µ T ) as given by Definition 3.1. Since we consider compact sets in finite dimensional spaces, both v and its gradient grad v can be simultaneously approximated uniformly by a polynomial up to an arbitrary precision. Therefore as f is bounded on X × U , there exists a polynomial v k of degree k such that sup

and finally, l ∈ polyIOCP(µ, µ T ) because was arbitrary (fixed).

All properties of Lagrangians in Definition 3.1 hold for the Lagrangians in Definition 4.4. For example, as stated in Remark 2, if l ∈ polyIOCP ,k (µ, µ T ), then (µ, µ T ) is close to optimal for OCP(l, µ 0 , T ). In particular, if we the moments of µ and µ T are available then the latter property can be checked numerically by solving a semi-definite program.

Integral discretization.

For practical numerical computation in the context of Definition 4.4, we still must be able to integrate polynomials with respect to µ and µ T . This is easy provided that we know the moments of µ and µ T . However, exact computation of such moments cane be complicated in practice, especially when µ is a superposition of trajectories. Usually, data sets from experiments consist of samples of trajectories which can be seen as realizations of a random sampling process.

5.3. Illustration on more complex problems. These simulations are taken from [START_REF] Pauwels | Inverse optimal control with polynomial optimization[END_REF]. We consider the following free terminal time direct problems:

Minimum exit time in dimension 2:.

The optimal law is u = x ||x||2 and the value function is v 0 (x) = 1 -||x|| 2 . Minimum exit norm in dimension 2:.

The optimal law is u = x and the value function is v 0 (x) = 1 -||x|| 2 2 . Minimum time Brockett integrator:.

(5.4)

Recall that the Brockett integrator of nonlinear systems control is also known (up to a change of coordinates) as the unicycle or Dubins system, one of the simplest instance of a non-holonomic system in robotics, see e.g. [START_REF] Devon | Kinematic and dynamic control of a wheeled mobile robot IEEE[END_REF] for the connection. The optimal law and value function are described in [START_REF] Prieur | Robust optimal stabilization of the Brockett integrator via a hybrid feedback[END_REF]. Complementary details are found in Appendix B of [START_REF] Pauwels | Linear conic optimization for inverse optimal control[END_REF]. Data generation:. We consider the following settings: A problem (5.2) with samples from B 2 ; A' problem (5.2) with samples from B 2 \ 1 2 B 2 ; B problem (5.3) with samples from B 2 ; C problem (5.4) with samples from B 3 \ 1 2 B 2 × R . In all cases we fix the degree of l to 4 and that of v to 10.

Results:. The results for the four problems are presented in Figure 5.2. For all problems, l is of degree 4. Therefore, l is to be found in a space of dimension 70 for problems A, A' and 126 for problem D. When the estimation error is close to 1, we estimate a Lagrangian l that is orthogonal to l 0 (in the monomial basis), and when it is close to 0, they are colinear. We also display the value of (iocp). We consider that the estimation is reasonable, when both estimation error and values are low.

• For all problems we are able to recover the true Lagrangian with good accuracy for some value of the regularization parameter λ. In the absence of regularization, we do not recover the true Lagrangian at all. This highlights the important role of 1 regularization which allows to bias the estimation toward sparse polynomials. The choice of λ in practical settings is subject to heuristics: numerical simulations or cross-validation which consists in keeping a portion of the input data as a validation set.

• For all four problems, when the estimation error is minimal, the value of is reasonably low, depending on how the value function can be approximated by a polynomial. For example, A' shows lower value because we avoid sampling database points close to the non-differentiable point of the true value function. In example D, the value function is known to be harder to approximate by polynomials and the value of is a bit larger. The estimation accuracy is still very reasonable.

Conclusion.

The main contribution of this paper is to propose a general framework to analyze the inverse problem of optimal control. The analysis is based on the weak formulation of direct optimal control problems using occupation measures, relaxed Hamilton-Jacobi-Bellman optimality conditions, and duality of infinitedimensional linear programs. The proposed formulation is powerful enough to ensure Appendix A. Proof of Proposition 4.7. We develop uniform finite sample bounds that hold with high probability for arbitrary probability distribution in the context of polynomial functions. These are in particular useful to derive bounds for the random process described by occupation measures as exposed in Section 4.3. The techniques used have become fairly standard in empirical process theory and statistical learning theory, see for example [START_REF] Bousquet | Introduction to statistical learning theory[END_REF] for a nice introduction.

In what follows, we consider a compact set Z ⊂ R p with non empty interior. For a polynomial z ∈ R p → q ∈ R d [z] of degree d, c(q) denotes its coefficients in the monomial basis (of size p+d d ). Similarly for a point z ∈ Z, v(z) denotes the p+d d dimensional vector representing the evaluation of the corresponding monomials at z such that q(z) = c(q) • v(z) with the dot denoting the inner product. We consider the following set of polynomials K d (Z) = {q ∈ R d [z], q ≥ 0 on Z, Z q = 1} where Z is a closed subset of Z with nonempty interior. We fix an arbitrary probability distribution P on Z. We denote by the linear functional on the space C(Z) such that

Similarly for a sample of size n, S n = {z 1 , . . . , z n } ∈ Z n , drawn iid from P , we denote by n the linear functional on the space C(Z) such that

for any function f continuous on Z. Lemma A.1. For any 0 ≤ δ ≤ 1 and any q ∈ K d (Z), it holds with probability 1 -δ that

where

are finite quantities that only depend on X and d.

Proof. The proof combines standard arguments from statistical learning which we describe here for completeness. In the sequel, given a probability distribution P on Z, we use the notation then, for all ≥ 0, when {z i } n i=1 is drawn iid from a probability distribution P on Z, we have

where the expectation is taken over the random sample. An equivalent formulation is that for 0 ≤ δ ≤ 1, with probability 1 -δ, it holds that

We consider the following quantity F (z 1 , . . . , z n ) := sup

n , p .