G Carlier
email: carlier@ceremade.dauphine.fr

A Lachapelle
email: lachapelle@ceremade.dauphine.fr

A planning problem combining optimal control and optimal transport

Keywords:

We consider some variants of the classical optimal transport where not only one optimizes over couplings between some variables x and y but also over some control variables governing the (possibly stochastic) evolutions of these variables with time. Such a situation is motivated by an assignment problem of tasks with workers whose characteristics can evolve with time (and be controlled). We distinguish between the coupled and decoupled case. The coupled case is a standard optimal transport with the value of some optimal control problem as cost. The decoupled case is more involved since it is nonlinear in the transport plan. We fully treat the linear quadratic case and show that even in this case existence may fail for the decoupled problem.

Introduction

In the classical optimal transport problem, a planner has to find some coupling or transport plan between two nonnegative measures of equal total mass (say 1) so as to minimize average cost. The Monge-Kantorovich problem is then the linear program consisting in finding a cost minimizing measure having prescribed marginals on variables, say x and y. Such problems have been very much studied in recent years in particular for their various applications in PDE's, geometry, probability theory (see the recent books of Villani [START_REF] Villani | Topics in optimal transportation[END_REF], [START_REF] Villani | Optimal transport. Old and new[END_REF])... Interestingly, transport/assigment problems have their modern roots in planning problems (optimally transporting coal from mines to steel factories for instance) and the pionneering works of Kantorovich in the USSR and Koopmans in the West were awarded by the Nobel prize in economics in 1975. Typical applications are optimal allocation of resources or assignment problems; for other economic applications, we refer to the recent work of Ekeland [START_REF] Ekeland | An optimal matching problem[END_REF], [START_REF] Ekeland | Existence, uniqueness and efficiency of equilibrium in hedonic markets with multidimensional types[END_REF] and Chiappori, McCann and Nesheim [START_REF] Chiappori | Hedonic price equilibria, stable matching, and optimal transport: equivalence, topology, and uniqueness[END_REF] relating matching problems and hedonic equilibria to optimal transport theory.

If, to fix ideas, we think of a firm having to assign N workers of different skills x to N tasks of different difficulties y, the assignment problem simply consists in assigning a task to each worker (for instance assigning most difficult tasks to most skilled workers). It is however reasonable to think that, in addition to designing such an assignment, the firm has the possibility to act directly on the workers' types using some (costly) control variable like internal training. The same of course applies to the task type y that may change with time if, for example, some investments are made to make the task easier or quicker to execute. This is where optimal control comes into play. The problems we consider in the present paper are precisely intended to deal with such situations where a planner or firm has at disposal as decision variables not only the initial assignment but also control variables that affect the (possibly stochastic) dynamics of the state variables x and y. More precisely, we consider the case where the planner chooses the initial coupling between x and y, initially distributed according to respective probabilities µ 0 and ν 0 , and a control variable governing the evolutions of the characteristics x and y. It is assumed that a pair (worker/task) that is initially formed at time 0 with initial characteristics (x, y) will remain paired or married during a fixed period although its characteristics may evolve with time, this captures the idea that for instance workers are assigned a certain task during some minimal period.

We will consider two different cases that we respectively call the coupled and the decoupled case and will compare them. In the coupled case, one optimizes over initial transport plans and control for each initial pair (x, y). The coupled problem thus amounts to solve a standard optimal transport where the cost function simply is the value function of some optimal control problem. In the decoupled case, we will assume that y is constant in time and the control variable for x has to be chosen independently of y (for instance the planner offers the same training to all the workers of same skill, independently of the task they will be assigned). This case is much more involved since it is not linear (and not even lower semicontinuous) in the transport plan. Under quite general assumptions it turns out that the coupled and decoupled problems have the same value but that the latter may fail to have solutions (a counter-example will be given in section 4).

The paper is organized as follows. In section 2, the coupled problem is introduced, an existence result is proved and some cases of uniqueness are given. In section 3, we study the decoupled problem which, contrary to the coupled case, gives rise to a nonlinear optimization problem over transport plans. Finally, section 4 is devoted to a linear quadratic deterministic case that we treat in details and for which we compare the coupled and the decoupled case, we also give some counterexamples.

The coupled case 2.1 The optimization problem

We consider a finite horizon framework on the period [0, T ]. The state equation for the joint dynamics of

Z t = (X t , Y t ) is dZ t = b(Z t , u t )dt + σ(Z t , u t )dW t , t ∈ [0, T ] (2.1)
where W t is a standard 2d-dimensional Brownian motion on some filtered probability space and (u t ) t∈[0,T ] is an admissible control process i.e. it is a progressively measurable process with values in a control space U , assumed for simplicity to be a compact metric space. We will denote by U the set of admissible controls. We assume that b ∈ C(R 2d ×U, R 2d ), σ ∈ C(R 2d ×U, M 2d ) (M 2d denoting the space of 2d × 2d matrices with real entries) and that there is some constant K 1 such that for every (z

1 , z 2 , u) ∈ R 2d × R 2d × U one has |b(z 1 , u) -b(z 2 , u)| + |σ(z 1 , u) -σ(z 2 , u)| ≤ K 1 |z 1 -z 2 |. (2.2)
The Lipschitz assumption (2.2) guarantees the existence and uniqueness of a strong solution to (2.1) for every initial condition Z 0 = (x, y) and every admissible control process u ∈ U. We will denote Z x,y,u t this solution in the sequel.

To each initial condition (x, y) and each admissible control process u ∈ U is associated a cost

E T 0 c(Z x,y,u s , u s )ds + φ(Z x,y,u T ) .
We assume that both c and φ are nonnegative and continuous and that there exist a constant K 2 and an exponent k ≥ 1 such that for every (z

1 , z 2 , u) ∈ R 2d × R 2d × U : |c(z 1 , u) -c(z 2 , u)| + |φ(z 1 ) -φ(z 2 )| ≤ K 2 |z 1 -z 2 |(1 + |z 1 | k-1 + |z 2 | k-1 ). (2.3)
The last ingredient in the problem is a pair of Borel probability measures µ 0 and ν 0 on R d giving respectively the distributions of the initial state for x and for y. The initial distribution γ for the joint process Z is chosen by the controller, it is a coupling between µ 0 and ν 0 and thus belongs to the set of transport plan Π(µ 0 , ν 0 ) consisting of all Borel probability measures γ on R d × R d having µ 0 and ν 0 as marginals. We shall finally assume the following on µ 0 and ν 0 :

R d |x| k dµ 0 (x) + R d |y| k dν 0 (y) < +∞.
(2.4)

The controller's optimization problem then reads as inf

γ∈Π(µ 0 ,ν 0 ) R d ×R d inf u∈U E T 0 c(Z x,y,u s , u s )ds + φ(Z x,y,u T ) dγ(x, y) (2.5)
which can simply be rewritten as the optimal transport problem:

inf γ∈Π(µ 0 ,ν 0 ) R d ×R d v(0, x, y)dγ(x, y). (2.6)
where v is the value function

v(t, z) := inf u∈U E T t c(Z z,u s , u s )ds+φ(Z z,u T ) | Z z,u t = z , ∀(t, z) ∈ [0, T ]×R 2d .
We won't discuss whether the previous optimal control admits a solution (for each fixed z). This issue relies on whether the value function is regular enough so that one can define an optimal feedback control, this is of course tightly related with the regularity for the corresponding Hamilton-Jacobi-Bellman equation. Instead, we will focus on the optimal transport formulation (2.6).

Existence

It is now easy to check that the infimum in (2.6) is attained, although the arguments are very standard in optimal transport (see [START_REF] Villani | Topics in optimal transportation[END_REF]), we give a proof for the sake of completeness.

Theorem 2.1. Under the assumptions of paragraph 2.1, the value of (2.6) is finite and attained by some optimal transport plan.

Proof. First, it follows from assumption (2.3) and standard arguments in control theory (see for instance [START_REF] Fleming | Controlled Markov Processes and Viscosity Solutions[END_REF], section IV.6) that v(0, .) has polynomial growth and is continuous. More precisely, there is a constant

K 3 such that 0 ≤ v(0, x, y) ≤ K 3 (1 + |x| k + |y| k ), ∀(x, y) ∈ R 2d (2.7)
and, for all (z

1 , z 2 ) ∈ R 2d × R 2d : |v(0, z 1 ) -v(0, z 2 )| ≤ K 3 |z 1 -z 2 |(1 + |z 1 | k-1 + |z 2 | k-1 ). (2.8)
It immediately follows from (2.7) and (2.4) that the value of (2.6) is finite. Since Π(µ 0 , ν 0 ) is tight and v(0, .) is continuous, for each n ∈ N * , there is some

γ n ∈ Π(µ 0 , ν 0 ) such that R 2d min(v(0, .), n)dγ n ≤ R 2d min(v(0, .), n)dγ, ∀γ ∈ Π(µ 0 , ν 0 ). (2.9)
Again, by tightness of (γ n ) n , we may assume, taking a subsequence if necessary, that γ n converges narrowly to some γ ∈ Π(µ 0 , ν 0 ). We thus have, for every M > 0 and every γ ∈ Π(µ 0 , ν 0 ):

R 2d min(v(0, .), M )dγ = lim n R 2d min(v(0, .), M )dγ n ≤ lim sup n R 2d min(v(0, .), n)dγ n ≤ lim sup n R 2d min(v(0, .), n)dγ = R 2d
v(0, .)dγ and then letting M → +∞ we conclude that γ solves (2.6).

Some cases of uniqueness

Our aim now is to exhibit some special cases where there is uniqueness of an optimal transport plan for (2.6) and this plan is given by a transport map, that is of the form γ = δ S(x) ⊗µ 0 for some measure preserving map S between µ 0 and ν 0 (which corresponds to a deterministic optimal coupling).

It is indeed well-known (see [10], [START_REF] Levin | Abstract cyclical monotonicity and Monge solutions for the general MongeKantorovich problem[END_REF] [3]) that if the transportation cost v(0, ., .) satisfies a certain structural condition called the generalized Spence-Mirrlees or twist condition and if µ 0 is regular enough (absolutely continuous for instance) then there is uniqueness of an optimal transportation plan that is in fact given by a transport map. We won't insist here on the twist condition and the regularity required on µ 0 in general. Let us indicate some particular cases instead:

• d = 1, v(0, ., .
) of class C 2 (say) and the cross-derivative ∂ xy v(0, x, y) is everywhere strictly positive, in this case, as soon as µ 0 is atomless, there is a unique transport plan that is given by the unique nondecreasing map (or monotone rearrangement) pushing µ 0 forward to ν 0 ,

• v(0, x, y) = c 0 (xy) with c 0 strictly convex and differentiable (see [START_REF] Gangbo | The Geometry of Optimal Transportation[END_REF]) and µ 0 absolutely continuous and compactly supported, then there is also a unique optimal transport plan and the latter is given by a transport map,

• of particular interest is the quadratic case c 0 (z) = |z| 2 which was first solved by Brenier [START_REF] Brenier | Polar factorization and monotone rearrangement of vectorvalued functions[END_REF], in this case, uniqueness of an optimal plan, once again given by a transport map, holds as soon as µ 0 and ν 0 have second moments and µ 0 does not give mass to sets of Hausdorff dimension less than d -1 (see McCann [START_REF] Mccann | Existence and uniqueness of monotone measurepreserving maps[END_REF]).

Uncontrolled case Let us first consider the uncontrolled case where b = b(z) is a Lipschitz vector field, σ = √ 2I 2d , the running cost c = c(z) and the terminal cost φ = φ(z) satisfy the same assumptions as in section 2.1. In this case, it is well-known that the function

v(t, z) := E T t c(Z s )ds + φ(Z T ) | Z t = z , (t, z) ∈ (0, T ) × R 2d
is the unique solution with polynomial growth of the backward linear parabolic equation

∂ t v + ∆v + b • ∇v = -c, v| t=T = φ. (2.10)
As seen in the previous paragraph, in the optimal transportation problem (2.6), the cost to take into account is v(0, .), the value at time 0 of the solution of (2.10). Let us treat in details the case b = 0 and φ = 0. In this case, (2.10) is simply the heat equation

∂ t v + ∆v = -c, v| t=T = 0, so that v(0, .) = T 0 (G t ⋆ c)(.)dt with G t (z) = 1 (4πt) d/2 exp - |z| 2 4t . (2.11)
Then the optimal transportation problem (2.6) reads as:

inf γ∈Π(µ 0 ,ν 0 ) R d ×R d T 0 (G t ⋆ c)(x, y)dt dγ(x, y).
In the particular case of the quadratic transportation cost c(x, y) = |x-y| 2 2 , a direct computation leads to:

v(0, x, y) = T |x -y| 2 2 + 1 .
Consequently, the infimum of (2.6) is T (W 2 (µ 0 , ν 0 ) + 1) where W 2 stands for the squared-2-Wasserstein distance:

W 2 (µ 0 , ν 0 ) := inf γ∈Π(µ 0 ,ν 0 ) R d ×R d |x -y| 2 dγ(x, y).
It is well-known that for the quadratic cost, if µ 0 does not give mass to small sets (i.e. Borel sets of Hausdorff dimension less than d -1) then there is a unique optimal transport plan that is in fact induced by a transport map S. Moreover S is characterized by S = ∇u with u convex (see [START_REF] Brenier | Polar factorization and monotone rearrangement of vectorvalued functions[END_REF], [START_REF] Gangbo | An elementary proof of the polar factorization of vectorvalued functions[END_REF], [START_REF] Mccann | Existence and uniqueness of monotone measurepreserving maps[END_REF], [START_REF] Gangbo | The Geometry of Optimal Transportation[END_REF]).

Linear control and zero transportation cost Let us now consider the linear drift b(Z t , u t ) = u t with constant diffusion σ = √ 2I 2d and the quadratic running cost c(z, u) = c(u) = u 2 2 . We also assume that φ ∈ C 2 (R 2d ) and satisfies the same assumptions as before. The value function:

v(t, z) := inf u∈U E T t c(u s )ds + φ(Z T ) | Z t = z , (t, z) ∈ (0, T ) × R 2d
is then the viscosity solution of the Hamilton-Jacobi-Bellman equation:

∂ t v + ∆v + |∇v| 2 2 = 0, v| t=T = φ.
(2.12)

Performing the usual Hopf-Cole transformation V (Tt, z) = e v(t,z)/2 , Φ := e φ/2 , transforms (2.12) into the heat equation for V :

∂ t V -∆V = 0, V | t=0 = Φ. (2.13)
As soon as Φ has polynomial growth, this gives the explicit expression

v(t, .) = 2 log(G T -t ⋆ Φ(.))
where G is again the heat kernel as in (2.11). In the unidimensional case, one may thus check whether v(0, ., .) satifies the cross-derivative criterion recalled above. Indeed, by a direct compuation, ∂ xy v(0, ., .) > 0 is equivalent to

(G T ⋆ ∂ xy Φ)(G T ⋆ Φ) > (G T ⋆ ∂ x Φ)(G T ⋆ ∂ y Φ).

Calculus of variations

Let us now consider the simple deterministic case where the variable y is constant with respect to time, where the state equation for X is simply Ẋ = u and where the running and terminal costs respectively take the form c(xy) + F (u), and φ(xy) where c, φ and F are bounded from below, smooth and strictly convex functions, with Λid ≥ D 2 F ≥ λid for some Λ > 0, λ > 0. Let us further assume for simplicity that µ 0 and ν 0 are compactly supported. In this case the value function is given by

v(0, x, y) = inf X∈H 1 , X(0)=x T 0 c(X(t) -y) + F ( Ẋ(t)) dt + φ(X(T ) -y).
We then rewrite the value as v(0, x, y) = w(xy) where

w(z) = inf X∈H 1 , X(0)=z T 0 c(X(t)) + F ( Ẋ(t)) dt + φ(X(T )).
It is now easy to check that w is a C 1 and strictly convex function. Then, it is well-known (see [START_REF] Gangbo | The Geometry of Optimal Transportation[END_REF], [START_REF] Gangbo | An elementary proof of the polar factorization of vectorvalued functions[END_REF], [START_REF] Carlier | Duality and existence for a class of mass transportation problems and economic applications[END_REF]) that the corresponding optimal transport problem inf γ∈Π(µ 0 ,ν 0 )

w(xy)dγ(x, y) possesses a unique solution that is given by an optimal transport map S as soon as µ 0 is absolutely continuous with respect to the Lebesgue measure (see [START_REF] Gangbo | The Geometry of Optimal Transportation[END_REF]). In this simple case again, we then have existence and uniqueness of an optimal plan that is further given by a transport map. We will again use this example in the decoupled variant of section 3.

Alternative formulations

We now describe formally an alternative formulation of (2.5). Given γ 0 ∈ Π(µ 0 , ν 0 ) as initial coupling of X and Y and an admissible control that we assume to be of Markovian type i.e. of the form u t = u t (X t , Y t ), let t ∈ [0, T ] and γ u t be the joint probability law of (X t , Y t ) (this is purely formal since the dependence of the Markovian control u on x, y may not be regular enough to define a flow to the SDE (2.1)):

R 2d ϕdγ u t := R 2d E(ϕ(Z x,y,u t ))dγ 0 (x, y), ∀ϕ ∈ C ∞ c (R 2d ).
Then γ u t (again formally) solves the Fokker-Planck equation (controlled by u):

∂ t γ - 1 2 tr σ(z, u t (z))σ(z, u t (z)) * D 2 γ + div(b(z, u t (z))γ) = 0, γ| t=0 = γ 0 .
(2.14) The total cost becomes linear with respect to γ u and simply reads as:

T 0 R 2 d c(z, u t (z))dγ u t (z) dt + R 2 d φ(z)dγ u T (z) (2.15)
which suggests as an alternative formulation of (2.5) to minimize (2.15) over initial conditions γ 0 ∈ Π(µ 0 , ν 0 ) and controls u, the dynamics of γ u being governed by the PDE (2.14).

Let us also mention that it is well-known that (2.6) admits the dual formulation: sup

R d ϕdµ 0 + R d ψdν 0 : ϕ(x) + ψ(y) ≤ v(0, x, y) .
(2.16)

The decoupled case

In the coupled problem considered previously, the situation is rather simple since the problem amounts first to compute the value function and then to solve the optimal transport problem where the cost is precisely this value function. However, this requires that the optimizer is able to design a control that depends on both state variables x and y. If one thinks of x as being workers' skills and y the difficulty of some task (constant in time to fix ideas), it may be reasonable to assume in certain cases that the control governing the evolution of x (training, education...) is the same for every worker of type x independently of the task y he will be assigned. In this case, one faces a somehow decoupled variant of the previous problem.

For the sake of simplicity, we will only consider in this section the deterministic case where in addition y is constant in time and Ẋ = u (calculus of variations). The running cost is given by some continuous Lagrangian (x, u, y) ∈ R 3d → L(x, u, y) that is convex with q-growth in the variable u for some q > 1, and the terminal cost is (x, y) ∈ R 2d → φ(x, y). We also assume for simplicity that µ 0 and ν 0 are compactly supported. In this framework, the coupled problem considered previously reads as (P) inf γ∈Π(µ 0 ,ν 0 )

I(γ) = R d ×R d v(x, y)dγ(x, y) where v(x, y) := inf X∈W 1,q , X(0)=x T 0 L(X(t), Ẋ(t), y)dt + φ(X(T ), y).
Of course, (P) possessses solutions as soon as v is continuous (or more generally lower semicontinuous). In the decoupled case, X is not allowed to depend on y but has instead to minimize the conditional average of the cost. More precisely, for γ ∈ Π(µ 0 , ν 0 ), let us denote the disintegration of γ with respect to its first marginal as γ = γ x ⊗ µ 0 and then define for every x:

v γ (x) := inf X∈W 1,q , X(0)=x R d T 0 L(X(t), Ẋ(t), y)dt + φ(X(T ), y) dγ x (y).
Then, the decoupled problem is:

(Q) inf γ∈Π(µ 0 ,ν 0 ) J(γ) = R d v γ (x)dµ 0 (x).
The novel feature of this decoupled problem is that, contrary to the optimal transport problem (P), it is nonlinear. In fact, since v γ (x) is concave in γ x (as an infimum of linear functionals of γ x ), (Q) is a concave minimization problem over Π(µ 0 , ν 0 ). Existence of a minimizer is a real issue here since the criterion in (Q) may not be weakly lower semicontinuous (and in fact even the measurability of x → v γ (x) is not obvious and a counterexample to lower semicontinuity will be given in the next section). Our aim now is not to study existence or nonexistence issues in general but simply to remark that (P) turns out to be a kind of relaxation of (Q) under quite weak assumptions: Proposition 3.1. If the value function v is continuous and µ 0 is atomless, then one has inf(Q) = min(P).

Moreover, γ solves (Q) if and only if it solves (P) and

v γ (x) = R d v(x, y)dγ x (y) µ 0 -a.e..
In particular, if S is an optimal transport for (P), i.e. δ S(x) ⊗ µ 0 solves (P), then it also solves (Q).

Proof. Since v γ (x) ≥ R d v(x, y)dγ x (y), one immediately deduces that inf(Q) ≥ min(P). Now, let us remark that if γ is given by a transport map that is γ = δ S(x) ⊗ µ 0 then v γ (x) = v(x, S(x)) for µ 0 -a.e. x and then in this case

one has R d v γ dµ 0 = R d v(x, S(x))dµ 0 (x) = R d ×R d vdγ.
Since µ 0 is atomless the set of transport plans induced by transport maps is weakly- * dense in Π(µ 0 , ν 0 ) (see [START_REF] Ambrosio | Lecture notes on optimal transport problems[END_REF]). In particular, since v is continuous, there is a minimizing sequence for (P) of the form γ n = δ Sn ⊗ µ 0 where S n ♯µ 0 = ν 0 for every n.

One then has

inf(Q) ≤ v γn dµ 0 = v(x, S n (x))dµ 0 (x)
and passing to the limit one gets inf(Q) ≤ min(P). Let γ ∈ Π(µ 0 , ν 0 ), γ solves (Q) if and only if

J(γ) = inf(Q) = min(P) ≥ I(γ)
that is γ solves (P) and

I(γ) = J(γ) i.e. v γ (x) = R d v(x, y)dγ x (y) µ 0 -a.e..
If we consider again the case

L(x, u, y) = c(x -y) + F (u), φ(x, y) = φ(x -y)
under the same convexity, growth and regularity assumptions as in the example of subsection 2.3, then we immediately deduce that if µ 0 is absolutely continuous with respect to the Lebesgue measure, then the corresponding decoupled problem admits a unique solution that is given by a transport map.

The linear-quadratic deterministic case

In this paragraph, we treat a special case where the state variable Y remains constant and the control on X is u = Ẋ. We will consider the case of a quadratic running cost and zero terminal cost:

c(x, y, u) = 1 2 |x -y| 2 + 1 2 |u| 2 , φ(x, y) = 0
and we will distinguish the coupled case (where u is allowed to depend on both x and y) and the decoupled case where u depends on the initial condition x but has to be the same for every y. Of course, in this quadratic framework, we will also assume that µ 0 and ν 0 have finite second-order moments.

The coupled case

In the coupled case, one first solves for every (x, y):

v(x, y) := inf

X∈H 1 (0,T,R d ), X(0)=x T 0 1 2 |X(t) -y| 2 + 1 2 | Ẋ(t)| 2 dt.
By standard arguments, this problem has unique solution given by:

X(t, x, y) := a(t)(x -y) + y, with a(t) = e t-2T + e -t 1 + e -2T , t ∈ [0, T ]. (4.1)
Replacing in v then leads to the optimal transportation problem:

1 2 T 0 (a 2 + ȧ2 ) inf γ∈Π(µ 0 ,ν 0 ) R d ×R d |x -y| 2 dγ(x, y) (4.2)
which is simply the optimal transportation problem with quadratic cost. Denoting by W 2 (µ 0 , ν 0 ) the squared-2-Wasserstein distance between µ 0 and ν 0 , we then have the following expression for the value of (4.2):

1 2 a 2 H 1 W 2 (µ 0 , ν 0 ) with a 2 H 1 = 1 -e -2T 1 + e -2T . (4.3)

The decoupled case

In the decoupled case, one has to minimize with respect to γ ∈ Π(µ 0 , ν 0 ) and (t, x) → X(t, x) (independent of y) such that X(0, x) = x, the quantity

R d ×R d T 0 1 2 |X(t, x) -y| 2 + 1 2 | Ẋ(t, x)| 2 dt dγ(x, y). (4.4)
By the disintegration theorem, every γ ∈ Π(µ 0 , ν 0 ) admits a disintegration with respect to its first marginal µ 0 i.e. can be written as γ = γ x ⊗ µ 0 where x → γ x is a Borel family of probability measures (naturally interpreted as conditional probability of y given x). We will denote by g γ (x) the average of γ x :

g γ (x) := R d ydγ x (y), for µ 0 -a.e. x ∈ R d .
Using the disintegration of γ with respect to µ 0 and the conditional moment g γ , it is convenient to rewrite the cost (4.4) as:

R d T 0 1 2 |X(t, x)| 2 -X(t, x) • g γ (x) + 1 2 | Ẋ(t, x)| 2 dt dµ 0 (x) + T 2 R d |y| 2 dν 0 (y).
Minimizing in X(., x) for each fixed x then yields

X(t, x) = a(t)(x -g γ (x)) + g γ (x) (4.5)
where a is as in (4.1). Replacing in (4.4), then yields

v γ (x) = 1 2 a 2 H 1 |x -g γ (x)| 2 + T 2 R d |y -g γ (x)| 2 dγ x (y). (4.6)
We thus find that the least cost, as a function of γ 0 ∈ Π(µ 0 , ν 0 ) only, is:

1 2 a 2 H 1 R d |x| 2 dµ 0 (x) + T 2 R d |y| 2 dν 0 (y) + 1 2 a 2 H 1 -T R d |g γ | 2 dµ 0 -a 2 H 1 R d x • g γ dµ 0 . (4.7)
Recalling that a 2 H 1 is given by (4.3), one easily gets a 2 H 1 -T < 0 for every T > 0 (which is consistent with the concavity property of the decoupled problem as noticed in section 3) and thus finding an optimal transport plan amounts to:

inf γ∈Π(µ 0 ,ν 0 ) J(γ) = -A(T ) R d |g γ (x)| 2 dµ 0 (x) -2 R d ×R d x • ydγ(x, y) (4.8)
where A(T ) is the positive constant:

A(T ) := T a 2 H 1 -1 = e -2T (T + 1) + T -1 1 -e -2T .
Of course when µ 0 is absolutely continuous, there is an optimal transport for the quadratic cost and then by proposition 3.1 it is also the unique solution of (4.8). The particular and simple form of the functional J in (4.8) enables us to deduce various counterexamples and some particular simple cases of non existence when µ 0 is atomless but may charge small sets.

The functional J is not lower semicontinuous: To see this, let us take µ 0 and ν 0 equal and uniform on the segment [-1, 1], γ n = δ Sn ⊗ µ 0 where S n is a sequence of µ 0 -measure preserving maps such that γ n weakly- * converges to the independent coupling µ 0 ⊗ µ 0 , then |g γ | 2 dµ 0 = 0 and then lim sup

n |g γn | 2 dµ 0 = lim sup S n 2 L 2 > |g γ | 2 dµ 0
since otherwise S n would converge to 0 in L 2 and then one would have γ = δ 0 ⊗ µ 0 .

A counterexample to existence: In view of the previous example, one should in general suspect non existence for the minimization of J over Π(µ 0 , ν 0 ). Let us consider the two-dimensional case where µ 0 is uniformly distributed on the segment [0, 1] × {0} (it thus charges a small set but still is atomless) and let ν 0 be uniform on the square [0, 1] 2 . Assume that there exists a minimizer γ of the decoupled problem. In this case, since µ 0 is atomless, using proposition 3.1, the explicit values of v(x, y) and v γ (x) given by (4.6), we would then have for µ 0 -a.e.

x v γ (x) = 1 2 a 2 H 1 |x -g γ (x)| 2 + T 2 R d |y -g γ (x)| 2 dγ x (y) = R d v(x, y)dγ x (y) = 1 2 a 2 H 1 R d |x -y| 2 dγ x (y) = 1 2 a 2 H 1 |x -g γ (x)| 2 + R d
|yg γ (x))| 2 dγ x (y) .

Which implies that γ x should be a Dirac mass for µ 0 -a.e. x, it would thus imply that γ is induced by a transport map S = (S 1 , S 2 ). At this stage it is not possible to conclude directly since there exist measure preserving maps between µ 0 and ν 0 , and even continuous ones such as Peano's function (see for instance [START_REF] Stone | Measure-preserving maps in General topology and its relations to modern analysis and algebra[END_REF]). To derive the desired contradiction, one has to use further proposition 3.1 to deduce that S is also an optimal transport map between µ 0 and ν 0 for the quadratic cost. Thus, by well-known arguments (see [START_REF] Gangbo | The Geometry of Optimal Transportation[END_REF], [START_REF] Villani | Topics in optimal transportation[END_REF], [START_REF] Villani | Optimal transport. Old and new[END_REF]) there exists a convex function (x, y) ∈ [0, 1] 2 → ψ(x, y) such that (S 1 (x, 0), S 2 (x, 0)) ∈ ∂ψ(x, 0) for (Lebesgue)-a.e. x ∈ [0, 1] which easily implies that S 1 (., 0) is nondecreasing and since S 1 (., 0) preserves the Lebesgue measure this yields S 1 (x, 0) = x. Now since S(x, 0) = (x, S 2 (x, 0)) is Borel and measure-preserving between µ 0 and ν 0 we should have for every Borel and bounded ϕ : [0, 1] 2 → R, so that taking ϕ(x, y) = (y -S 2 (x, 0)) 2 yields the desired contradiction.

The discrete case, different values for the coupled and decoupled problems: Let us now consider the case where both µ 0 and ν 0 have finite support, in this case (4.8) is in fact a finite-dimensional problem and it is easy to check that it possesses solutions. However, because of atoms, proposition 3.1 does not apply and it is easy to find a case where the values are different. To see this, let us choose d = 1, µ 0 = δ 0 and ν 0 = 1 2 δ -1 + 1 2 δ 1 . The set of transport plans reduces to the single element γ = 1 2 δ (0,1) + 1 2 δ (0,-1) . Obviously, the transport plan γ is not induced by any transport map. In this example, the solution of both the coupled and the decoupled problem is the transport plan γ. A simple computation leads to find different values for the two problems. More precisely, for a fixed horizon T > 0, one gets the difference: inf(Q) -inf(P) = 1 2 T -||a|| 2 H 1 > 0. Note that in this very simple case, the difference between the decoupled and coupled problem values is increasing in the horizon T .


J Masbou

G Lamanna

S Rosier-Lees

H E S S Collaboration

HESS-II expected performance in the tens of GeV

Keywords: Event reconstruction, multivariate analysis, energy reconstruction

By the end of 2009 the four-telescopes-HESS-system will be upgraded to HESS-II with a new telescope with a 600 m 2 mirror area and a very high resolution camera. HESS-II will allow to lower the energy threshold from 100 GeV to about 30 GeV and enhance the HESS sensitivity.

Therefore AGNs at higher redshift could be detected and searches for new classes of very high energy gamma-ray emitters (pulsars, microquasars, GRB, and dark matter candidates) will also be possible.

The evaluation of the instrument performance is presented in term of sensitivity, energy and angular resolutions, based on Monte Carlo simulation, using a multivariate analysis.

I. INTRODUCTION

By the end of 2009 a new very large telescope will be included in the present HESS-I configuration [START_REF] Hofmann | Status of the H.E.S.S. Project[END_REF]. The data taking configuration will consist of the four initial telescopes (called Cherenkov Telescope 1, CT2, CT3 and CT4) added to the new telescope (CT5) [START_REF] Punch | S. for higher sensitivity and lower energy[END_REF]. Since this new telescope will open a new window in the energy range from 30 to 100 GeV, it is important to develop a powerful analysis in order to optimize the gamma-ray detection and reconstruction. The new analysis presented here, includes only the fifth telescope. The selection of the gamma-ray events as well as the hadron rejection are described. A new energy reconstruction method is also reported.

II. ENERGY DOMAINS

The overall trigger system of the HESS-II configuration will operate in various modes, allowing to trigger on two classes of events :

1) Mono mode : Only the largest telescope CT5 has triggered, dedicated for the lowest energies 2) Stereo mode : As the current HESS system with additional information from the central telescope (CT5). For the evaluation of the low energy domain performance, spectrum of gamma-ray showers have been simulated from 5 GeV to 100 TeV as well as 18 discrete energy points from 5 GeV to 12.5 TeV, at the zenith, assuming an optimal optical efficiency. Proton showers were also simulated from 30 GeV to 100 TeV, in the same conditions. The gamma-ray source is simulated on the optical axis, thus the projection of the source is located in the center of the camera; proton showers are simulated isotropically (no point source). All the simulations are carried out on a large radius of 550 m from the center of the camera and a simulated night sky background of 100 MHz has been used. In what follows the performance of the standalone CT5 telescope is presented here without considering complementary information, if any, from the other telescopes.

III. DISCRIMINANT VARIABLES

The proton rejection is a key issue, especially at lower energy, and a multivariate analysis is proposed to this end. The usual Hillas parameters [START_REF] Hillas | Cerenkov light images of EAS produced by primary gamma[END_REF] are used. In addition a new variable Δφ is included which takes into account the background morphology as detailed below. A preselection is applied for the events at the edge of the camera, their image being not fully contained: nominal distance 1 < 1 • . The raw image in the telescope is first cleaned; different cleanings of the image are available, depending on the applied threshold. A tight cleaning (5/10) 2 with a high threshold allow to reconstruct the core of the shower, a second one (1/3) with a low threshold is sensitive to the sub-structures. Assuming the shower image is an ellipse, the six Hillas parameters for the core 3 of the shower are computed. Δφ is the difference between the orientations of the two ellipses using the two different cleaning. Therefore Δφ is dependent on the hadron topology and will not change for the gamma signal while it will for hadrons. Additional checks have been performed concerning the night sky background dependence assuming different frequencies ranging from 30 to 300 MHz, as depicted in figure 1.

IV. MULTIVARIATE ANALYSIS FOR GAMMA ENERGY RECONSTRUCTION For a stereoscopic system, the charge and the impact parameter suffice to determine the energy. It is not the case when considering only one telescope. As it is well known the Hillas variables are sensitive to the energy, the fixed energy points have been used in what follows. 1 The nominal distance is defined as the distance between the center of the camera and the barycenter of the reconstructed shower image 2 Cleaning 5/10 corresponds to a charge threshold of 10 photoelectron (p.e.) and 5 p.e. to one neighbour. 3 length, width, skewness, kurtosis, total charge and nominal distance. A new method is proposed to measure the energy in order to reduce the bias. The neural network approach uses the module MultiLayerPerceptron developed in the ROOT framework (TMVA) [START_REF]Toolkit for Multivariate Data Analysis with ROOT[END_REF] in order to determine the most probable energy for each event.

1) Training step : The signal events are Monte Carlo gamma-ray in a given energy band, the other energy bands are considered as the background.

For two different trainings, there is an overlap to insure the continuity of the energy reconstruction. 2) Analysis step : For each event the neural network (MLP) allows to determine the probability of each energy band (corresponding to the 12 previous trainings). The reconstructed energy is found to be the most probable value (i.e. the peak) as depicted in figure 2. V. HADRON REJECTION Thus for each event, the reconstructed energy is computed, applying the previous algorithm developed for gamma-ray events. The gamma/hadron separation is performed after using the reconstructed energy and different discriminant variables previously defined (Hillas parameters and Δφ). A similar neural network frame is used for the selection and the hadron discrimination.

MC Energy [GeV]

1) Training step : The signal events are Monte Carlo gamma-ray in the same energy bands, the background is the Monte Carlo proton spectrum. 2) Analysis step : The reconstructed energy is used to select the appropriate neural network. The final optimization corresponds to a global background rejection of 95%.

VI. PERFORMANCE For the events passing the selections (total background rejection of 95%), the expected energy and angular resolutions, and the effective area for CT5 standalone are estimated. The effective area after selections varies from 40,000 m 2 at 30 GeV to 80,000 m 2 at 100 GeV lowering the energy threshold to around 30 GeV, roughly one order of magnitude lower than HESS phase I (Fig. 3).

The new energy measurement method seems to be linear as depicted in figure 4. The difference between the energy simulated and reconstructed (so called bias) is flat in the [30 -300] GeV range and lower than 5%. Thanks to the linearity, the energy resolution can now be estimated: it decreases from 40 % to 20 % for energies between 30 and 500 GeV as shown in figure 5.

The angular resolution, defined as the 68% containment radius, varies from 0.25 to 0.15 deg for 30 to 500 GeV energy range, as shown in figure 6.

VII. CONCLUSION

This new energy reconstruction method shows very promising results. In particular the energy measurement is linear within 5% in the [30 -300] GeV energy range. The effective area in the CT5 standalone mode after selections varies from 40,000 m 2 at 30 GeV to 80,000 m 2 at 100 GeV lowering the energy threshold to around 30 GeV. The energy resolution varies from 40 to 20% as a function of the energy. The angular resolution is also improved, better than 0.25 deg. Nevertheless, there is still some room for improvement such as exploiting the timing information, combining the data of all the telescopes or including the future HESS-II trigger.

Fig. 1 :Fig. 2 :

 12 Fig. 1: Distribution of Δφ assuming different night sky background from 30 to 300 MHz.

Fig. 3 :

 3 Fig.3: CT5 standalone: acceptance after selection (continuous blue line) compared to the trigger acceptance (upper continuous green line). Superimposed: HESS-I acceptance (dashed line).

Fig. 4 :Fig. 6 :

 46 Fig. 4: CT5 standalone: Energy linearity measurement.


Stéphane Weiss
email: weiss@loria.fr

Pascal Urso
email: urso@loria.fr

Pascal Molli
email: molli@loria.fr

Compensation in Collaborative Editing

Undo/Redo has been recognized as an important feature of collaborative editing systems. The undo mechanism allows any user to undo any edit operation at any time. Preserving consistency of shared data with the undo feature is a complex issue. Several algorithms exist for managing Undo in collaborative editing, but they are still open issues concerning correctness or intrinsic limitations of these algorithms. Therefore, in this paper we present an novel undo approach in the context of Operational Transformation (OT) framework. We prove the correctness of our solution, we leverage some limitations of current Undo approaches and we indicate how it can be combined with existing OT integration algorithms to provide a complete usable system.

INTRODUCTION

Undo/Redo has been recognized as an important feature of collaborative editing systems [START_REF] Abowd | Giving undo attention[END_REF][START_REF] Choudhary | A general multi-user undo/redo model[END_REF][START_REF] Sun | Achieving convergence, causality preservation, and intention preservation in real-time cooperative editing systems[END_REF]. In collaborative systems, the most general model of undo mechanism allows any user to undo any edit operation at any time. Preserving consistency of shared data with the undo feature is a complex issue. Several Undo algorithms have been published within the Operational Transformation (OT) [START_REF] Sun | Achieving convergence, causality preservation, and intention preservation in real-time cooperative editing systems[END_REF][START_REF] Ellis | Concurrency control in groupware systems[END_REF] framework.

Ressel proposed the first undo algorithm with Adopted [START_REF] Ressel | Reducing the problems of group undo[END_REF]. However, this undo algorithm does not allow the user to undo any operation. The user can only undo his proper operations in the inverse order of their generation.

Next, the ANYUNDO approach [START_REF] Sun | Consistency maintenance in real-time collaborative graphics editing systems[END_REF] provides a mechanism to undo any operations. However it requires the verification of the properties TP1, TP2 and IP1. Unfortunately, transformation functions satisfying these properties have never been published. The undo algorithm from [START_REF] Ferrié | Concurrent undo operations in collaborative environments using operational transformation[END_REF] suffers from the same limitations as the ANYUNDO approach.

The COT approach [START_REF] Sun | Operation Context and Context-based Operational Transformation[END_REF] has been designed to simplify the ANYUNDO approach. The major changes presented in COT include, along the integration algorithm, the usage of context vectors and the breaking of property TP2. In this paper, we demonstrate that property TP2 is not broken in the COT approach. We also point out that context vectors grow linearly with the number of undo performed during the editing session. Thus, context vectors limit the number of operation that we can undo. Consequently, the number of operation a user can undo depend on the number of participants and the network's bandwidth.

In this paper, we present the compensation approach. This approach does not require context vectors and property IP1. Moreover, compensation can be used with all integration algorithms such as adOPTed, SOCT2, GOTO and COT. It makes the undo feature available even for algorithms with no native undo support such as SOCT2 and SOCT4.

As a proof of concept, we applied the compensation framework to the Tombstone Transformation Functions (TTF) [START_REF] Oster | Tombstone transformation functions for ensuring consistency in collaborative editing systems[END_REF]. The correction of the obtained model is proved formally by the automated proof environment VOTE [START_REF] Imine | Vote: Group editors analyzing tool: System description[END_REF].

Based on this approach, we build the Graveyard realtime editor prototype. Graveyard combines a SOCT2 algorithm that does not provide a native undo feature, with TTF transformation functions extended with compensation operations.

COMPENSATION IN THE OT APPROACH

Existing undo methods for OT enforce the document to return to a previous state after the execution of an undo operation. The compensation approach does not require the system returns to a previous state but only to a state semantically equal to a previous state. However, this semantically equal state can also be syntactically equal to the previous state.

We call C(op) the operation which compensates op. The execution of the operation C(op) undoes, from a semantic point of view, the effect of the operation op. C(op) is not necessary the inverse operation of op but it is defined in order to compensate op just after the execution of op. We need to define each operation C(op) in such a way that S • op • C(op) is a state where the effect of op has been semantically undone.

However, the operation C(op) does not take account of the effect of operations executed between the execution of op and it's compensation by the user. We need to compute an operation C(op) ′ which realize the effect of C(op) on the current state. To compute the operation C(op) ′ , we need the following algorithm (also illustrated by the Figure 1).

1. First, we determine the operation C(op) which compensates op. C(op) should have compensated op if op was the last operation executed on the current site.

2. We use the forward transformation functions in order to transform C(op) with all operations which have already been executed on this site. The resulting operation is called C(op) ′ . C(op) ′ is defined on the current state and will be send to other sites where it will be integrated as any other new operation.

This algorithm is known as the naive algorithm for undo [START_REF] Ferrié | Concurrent undo operations in collaborative environments using operational transformation[END_REF].

site 1

C(op) op op 1 op 2 C(op) ′ = T (C(op), op 1 • op 2 )
Figure 1. Algorithm of the compensation

CORRECTNESS OF THE COMPENSATION APPROACH

The compensation approach does not depend of the integration algorithm. It is defined at the transformation functions level. In order to be correct, the transformation functions defined for regular and compensation operations have to ensure the following properties.

• The properties T P 1 and T P 2 are defined to ensure consistency. Depending on the OT integration algorithm, both properties are required or only T P 1. Indeed, in our approach, integration algorithms make no distinction between compensating operations and other operations. Then, compensating operations have also to satisfy the correctness properties.

• The property T P C ensures the respect of the compensation. This property ensures that the compensation effect will always be the same, even if the operation compensated is not the last executed one. This condition is similar to the conditions C4 [START_REF] Ferrié | Concurrent undo operations in collaborative environments using operational transformation[END_REF] and IP 3 [START_REF] Sun | Operation Context and Context-based Operational Transformation[END_REF][START_REF] Sun | Undo any operation at any time in group editors[END_REF].

T P C : T (C(op), T (seq, op)) = C(T (op, seq))
Figure 2 explains the T P C property. Two sites make concurrent operations. Site1 generates op while site 2 generates a sequence of operations seq. Both sites receive remote operations, transform and integrate them. Now, they are on the same state. Consequently, if they want to compensate the same operation on the same state, they must obviously generate the same operation. Site1 generate C(op) and transform it through following operations T (seq, op). Site2 compensate the last received operation which is T (op, seq). These two compensating operations are defined on the same state, they compensate the same operation, so they must be the same. The verification of this property ensures that whenever an operation is compensated, the compensation effect remains the same.

So, there are three -or two -properties to verify in order to ensure a correct OT system with compensation. Due to their conciseness, these properties are theoretically easy to prove. However, one of the particularity of the OT approach is the huge numbers of cases to check.

In such conditions, a hand proof is error-prone, and many transformation functions supposed hand-proven finally revealed themselves false (all counter examples can be found in [START_REF] Oster | Tombstone transformation functions for ensuring consistency in collaborative editing systems[END_REF]). On another hand, each of the cases to check can be easily handle by an automated formal theorem prover. Consequently, we choose to use the proof environment VOTE [START_REF] Imine | Vote: Group editors analyzing tool: System description[END_REF] based on the theorem prover Spike [START_REF] Stratulat | A general framework to build contextual cover set induction provers[END_REF] which generates all cases and ensures the verification of all required properties.

COMPENSATION IN THE TTF APPROACH

The TTF approach is divided in two parts: the model and the transformation functions. A detailed explanation of the TTF approach and its correctness can be found in [START_REF] Oster | Tombstone transformation functions for ensuring consistency in collaborative editing systems[END_REF].

In the TTF approach, deleted characters are kept as tombstones. The document's view only shows visible characters and tombstones are hidden. Consequently, the model differs from the view. Assume that a document is in a state "abcd". Now, a user deletes the character 'b'. In the TTF model, the character is replaced by a tombstone (i.e., the character with a visibility flag set to false). The view differs from the model as the view only contains "acd" while the model contains "ab /cd". Since tombstones are necessary to achieve consistency, they cannot be removed and thus, the operation "Ins" is not inversible.

As the TTF approach is defined for two operations: "Ins(p, c, sid)" and "Del(p, sid)", we need to find two operations to compensate them. The compensating effect required is to obtain that visible characters are the same.

Fortunately, it is obvious that there is no difference between compensating an insertion and deleting a character. So we can use the operation "Del" for compensating "Ins" since it has the desired effect which is removing the character in the user's view.

However, compensating the operation "Del" with the operation "Ins" inserts a new character in spite of the existing tombstone for this character. Thus, we define a new operation "Undel(position, sid)" to compensate the operation "Del". This operation replaces the tombstone of a deleted character by the original character. This operation "Undel" allows to build transformation functions which satisfy T P C and to reuse tombstones.

The function C(op) links normal operations to compen- Finally, we can write the transformation functions for all operations. The definition of the transformation functions for the operations "Ins" and "Del" are the same as presented in [START_REF] Oster | Tombstone transformation functions for ensuring consistency in collaborative editing systems[END_REF].

Site1 Site2 C(op) A A op o o C C Seq s s f f f f f f f f f f f f f f f f f f Seq ′ = T (Seq, op) op ′ = T (op, Seq) G G C(op ′ ) r r C(op) ′ = T (C(op), Seq ′ ) C(op ′ ) = C(T (op, Seq)

Correction of the approach

In the TTF approach, the transformation functions are written in order to satisfy the property T P 2. The property T P 1 is defined by a state equality, then we have to define the effect of the operation "Undel" on the state. If the effect of "Undel" is simply to make the character visible, the property T P 1 is violated (see Figure 3) To ensure T P 1, we define the effect of "Undel" replacing the visibility flag associated to characters by a visibility level. This visibility level is an integer. Initially, a character inserted has a visibility level of 1. Each time an operation deletes this character, its visibility level is decreased. Each time an operation undeletes this character, we increase its visibility level. A character is said Site1 "ab" Site2 "ab" Del(0, sid 1 )

8 8 x x x x x x x x x x x
Del(0, sid 1 )

Ð Ð Ñ Ñ Ñ Ñ Ñ Ñ Ñ Ñ Ñ Ñ Ñ Ñ Ñ Ñ Ñ Ñ Ñ Ñ "a /b" "a /b"
U ndel(0, sid 1 )

8 8 x x x x x x x x x x x Del(0, sid 1 )
"ab" "a /b" Del(0, sid 1 ) U ndel(0, sid 1 ) "a /b" "ab" "visible" and appears in the document's view if its visibility level is at least 1. Similarly, a character is said "invisible" and does not appear in the document's view if its visibility level is less than 1 (see figure 4).

Site1 "a 1 b 1 " Site2 "a 1 b 1 "
Del(0, sid 1 )

8 8 x x x x x x x x x x x Del(0, sid 2 ) Ð Ð Ñ Ñ Ñ Ñ Ñ Ñ Ñ Ñ Ñ Ñ Ñ Ñ Ñ Ñ Ñ Ñ Ñ Ñ "a / 0 b ′′ 1 "a / 0 b ′′ 1 U ndel(0, sid 1 ) 8 8 x x x x x x x x x x x Del(0, sid 1 ) "a 1 b ′′ 1 "a / -1 b ′′ 1 Del(0, sid 2 ) U ndel(0, sid 1 ) "a / 0 b ′′ 1 "a / 0 b ′′ 1 Figure 4. Visibility level
Using the proof environment VOTE [START_REF] Imine | Vote: Group editors analyzing tool: System description[END_REF], we have proven that our transformation functions satisfy the properties T P 1, T P 2 and T P C . The system specification given to the theorem prover Spike can be reviewed and tested at the following url : http://potiron.loria. fr/projects/graveyard.

Implementation

In order to validate our approach, we have built the Graveyard prototype. Graveyard is a real-time collaborative text editor. It relies on the SOCT2 algorithm for integrating concurrent operations. SOCT2 does not provide natively undo capabilities. We used the TTF transformation functions with related compensation operations to obtain a real-time collaborative with an undo feature. The general architecture of graveyard is described in figure 5.

For this implementation, we used SOCT2 but we can replace SOCT2 by SOCT4, adOPTed or COT and obtain the same result. Indeed, our framework only introduces the compensation algorithm. To compensate an operation op, we transform C(op) with all operations which have been executed after op. Fortunately, the main goal of every integration algorithms is to transform an operation against a set of concurrent operations. Thus, every OT system already contains such a feature : Translate Request for adOPTed, Transpose Forward for SOCT4 or COT-DO for COT. Finally, the compensation approach can be used any integration algorithm.

RELATED WORK

In [START_REF] Ressel | Reducing the problems of group undo[END_REF], the authors present an undo specific to the adOPTed algorithm. Unfortunately, this solution cannot allow to undo any operation at anytime.

The ANYUNDO algorithm [START_REF] Sun | Consistency maintenance in real-time collaborative graphics editing systems[END_REF] is associated with the GOTO integration algorithm. The GOTO-ANYUNDO approach needs transformation functions which satisfy three properties T P 1, T P 2 and IP 1. The property IP 1 illustrates the neutrality of do-undo pairs toward the document state. Unfortunately, transformation func-tions satisfying all these properties have never been published. For instance, TTF functions satisfy T P 1 and T P 2 but not IP 1.

In [START_REF] Ferrié | Concurrent undo operations in collaborative environments using operational transformation[END_REF], the authors define two properties C3 and C4 which are similar to IP 2 and IP 3. To ensure the verification of these two properties, the authors introduce a specific operation "undo(op)". This approach defines generic transformation functions for this operation "undo(op)" using the proposed transformation functions.

Unfortunately, as the ANYUNDO, a property similar to IP 1 is required. Therefore, this approach cannot be instantiate as it required properties which have never been published.

In COT [START_REF] Sun | Operation Context and Context-based Operational Transformation[END_REF], the authors propose an approach which requires only T P 1. To ensure consistency, concurrent operations are ordered before transforming remote operations. Unfortunately, this ordering is not sufficient as shown by Figure 6. Let's assume that concurrent operations order is op 1 , op 2 and op 3 . When site3 receives op 1 , due to concurrent operations order, op 1 is transformed with op 2 • op ′ 3 with op ′ 3 = T (op 3 , op 2 ). Similarly, site4 transforms op 2 with the sequence op 1 • op ′′ As the property TP2 is not broken, we need transformation function satisfying TP1 and TP2. Therefore, we apply the TTF functions with the COT algorithm (Figure 7) Moreover, the COT approach introduces the notion of context vector. A context vector is composed by a classical state vector and an extra component designed for undo. This extra component stores all undone operations. Each of them is associated with an integer representing how many times this operation has been undone. This context vector is associated to each operation. Its size grows linearly with the number of undone operation. This means that the upload bandwidth will be eventually saturated by the growing context vectors.

To illustrate this issue, let's assume that 10 users are editing the same text document using the COT algorithm. They are connected through Internet using an ADSL connexion. The upload limit for a such connexion is usually 1M b/s. We assume that a user generates 4 operations per second. If 32bits integer are used, the bandwidth is saturated by context vectors after 4091 undo 2 . Consequently, each user can only undo around 400 operations and, thus, affects the same number of characters.

On another hand, the TTF approach, which is independent of the compensation approach, the document model keeps the tombstones and is thus also unbounded. However, the operations sent on the network have a

does not affect the result. 2 We only count context vectors without any encapsulations such as IP,TCP,... fixed size and the bandwidth consumption is limited.

Finally, due to the space complexity of context vector, we consider that the COT approach only supports few users connected on an high bandwidth network.

CONCLUSIONS

Many approaches aimed to provide an undo mechanism which allows a user to undo any operations at any time.

Unfortunately, most of them are not instantiable and the other one generates heavy traffics which are not compatible with real-time constraints.

In this paper, we introduced our compensation mechanism. The compensation approach is more generic than existing undo approach: we can apply compensation to all transformation functions even if some operations have no inverse. An important feature of our approach is that the resulting transformation functions remain generic towards integration algorithms. Consequently, we can apply these functions with COT, SOCT2, SOCT4, GOTO and adOPTed. We have a complete solution to build text editors with undo capabilities. The compensation approach proposed in this paper has been implemented in the Graveyard collaborative text editor based on the tombstone transformation approach.

In future works, we will implement more existing integration algorithms in the Graveyard prototype. Thus, we will be able to benchmark the different integration algorithms.

Figure 2 .

 2 Figure 2. Respect of the compensation effect

Figure 3 .

 3 Figure 3. Violation of the property T P 1

Figure 5 .

 5 Figure 5. Graveyard architecture

3 2 =Figure 6 .

 26 Figure 6. COT Counter-example.

op 2 = 3 =Figure 7 .

 237 Figure 7. Using TTF functions with COT.

properties IP2 and IP3, but the COT approach aims to enforces them. In figure 7, when site1 receives op 2 , it has to transform this operation with op 1 • op 3 . As this sequence is a do-undo-pair, the COT algorithm skips this sequence and it does not transform op 2 . Using TTF with COT does not ensure consistency in case of undo. Thus, the COT approach fails to enforce IP2 in some cases. And finally, none of the published transformation functions can be directly used with COT to build a text editor with undo capability. site 1 "abc" site 2 "abc"


Sébastien Lengagne

Nacim Ramdani

Philippe Fraisse

Planning and Fast Re-Planning of Safe Motions for Humanoid Robots : Application to a Kicking Motion

Optimal motions are usually used as joint reference trajectories for repetitive or complex motions. In the case of soccer robots, the kicking motion is usually a benchmark motion, computed off-line, without taking into account the current position of the robot or the direction of the goal. Moreover, robots must react quickly to any situation, even if not expected, and cannot spend time to generate a new optimal motion by the classical way. Therefore, we propose a new method for fast motion re-planning based on an off-line computation of a feasible sub-set of the motion parameters, using Interval Analysis.

INTRODUCTION

Planning optimal motions for humanoid robots is dedicated to complex or repetitive tasks. Humanoid robots are complex systems in which neither geometric structure nor dynamic model are simple. This complexity limits the reactive capabilities of computation for motion planning. Indeed, a set of constraints such as balance, maximal joint torque velocity or position has to be included into the optimization process.

These constraints can be non linear and require a large computation time. Consequently, optimal motions are often generated off-line and used as joint reference trajectories. Motion planning includes as well, the problem of digital actors' locomotion [START_REF] Yoshida | Humanoid motion planning for dynamic tasks[END_REF], kick motion generation on HRP-2 robot [START_REF] Miossec | Development of a software for motion optimization of robots -application to the kick motion of the hrp-2 robot[END_REF], computing a manipulator robot's trajectory [START_REF] Piazzi | Global minimum-jerk trajectory planning of robot manipulators[END_REF] or smoothing pre-calculated motions [START_REF] Suleiman | On humanoid motion optimization[END_REF]. In previous works, we presented a new method for planning safe motions [START_REF] Lengagne | Guaranteed computation of constraints for safe path planning[END_REF], [START_REF]A new method for generating safe motions for humanoid robots[END_REF], [START_REF]Safe motion planning computation for databasing balanced movement of humanoid robots[END_REF] which uses the Interval Analysis to compute the constraint functions over time-interval, whereas classical methods which compute them over a time grid without any information about the constraint validity between two points of the grid.

Nowadays, motion planning aims to be fast enough to react to unexpected events. For instance, homotopic paths allow to modify a generated trajectory, to address collision avoidance for mobile arm manipulator [START_REF] Brock | Real-time re-planning in high-dimensional configuration spacesusing sets of homotopic paths[END_REF].

In the case of humanoid robots, the computation of the inequality constraint functions is too time-consuming, and motion planning can last from a few minutes to several hours. However, all the motions planned, or re-planned must satisfy a set of inequality constraint functions (balance, maximum torque or angle joint).

The re-planning process starts from an optimal motion parameter set, and computes a feasible sub-set that is an inner approximation of the feasible set of the inequality constraints. Then, it has to find a solution, in this subset, which allows to overcome unpredicted situations. By the way, the re-planning process consider the inequality constraints without computing them online.

In the case of soccer robots [START_REF] Behnke | See, walk, and kick: Humanoid robots start to play soccer[END_REF], the most important motion is the kicking motion, since it allows to goal. Usually, kicking motion is computed off line [START_REF] Carpin | The challenge of motion planning for human robots playing soccer[END_REF], hence does not take into account the current position of the robot or the direction of the goal. Nevertheless these motions allow the robot to react quickly to the situation, even if the kicking can lead to an accurate trajectory of the ball. In this paper, we propose a method to make the kicking motion more accurate by an off-line planning and a fast re-planning process.

In Section I, we briefly show how to generate optimal motions. Then, Section II presents our algorithm for computing a feasible box around the optimal motion. We apply this method to a kicking motion for the HOAP-3 Humanoid robot in Section III

I. SAFE MOTION PLANNING

A. Modeling 1) B-splines parameterization: Motion planning problem is an infinite programming problem which can be transformed into a Semi-Infinite Problem (SIP) by a joint trajectory parameterization [START_REF]A new method for generating safe motions for humanoid robots[END_REF].

We choose to compute the joint trajectories thanks to Bsplines functions [START_REF] De Boor | A Pratical Guide to Splines[END_REF]. Thus, we define a motion via the parameter vector X = [T, p 1,1 , p 1,2 , . . . , p 6,5 ] where T is the motion duration and p i the coefficients of the weighted sum of B-splines functions. The joint trajectory q i (t) is computed as follows: (./#4428#./$+393/:$1.;$,<:/#=: >7/39#0$??@?AB$%&&'$,/C$D3E2:B$F,G

q i (t) = N s ∑ j=0 p i, j × B j (t) (1)
The joint velocity and acceleration are obtained by differentiating the equation [START_REF] Yoshida | Humanoid motion planning for dynamic tasks[END_REF]. In this paper, we consider 5 Bsplines functions represented in Figure 1. Those functions allow to get motion with initial and final velocity and acceleration equal to zero.

2) Inverse Dynamic Model: We model the humanoid robot as an arborescent chain with n degrees of freedom. We assume that the right foot is the reference body, and add some constraint during the motion planning to make sure that is motionless. We use HuManS software [START_REF] Wieber | The humans toolbox, a homogeneous framework for motion capture, analysis and simulation[END_REF] to generate the models in the form of C-files functions. Thanks to the definition of Denavit-Hartenberg parameters [START_REF] Khalil | Modeling, Identification & Control of Robots[END_REF], center of mass and inertia values, this software generates C-functions that compute the inertial effect matrix M(q) and the Nonlinear effects vector due to gravity and Coriolis H(q, q).

Starting from the vectors {q(t) ∈ R n+6 , q(t) ∈ R n+6 , q(t) ∈ R n+6 } which contain the angle value of the n joints and the position of the reference body, we compute the joint torques, Γ(t) ∈ R n , and the effort applied by the reference body (the right foot) to the ground

F(t) ∈ R 6 , Γ(t) F(t) = M(q(t)) q(t) + H(q(t), q(t)) (2)
The HOAP-3 Humanoid Robot, represented in Figure 2, has 21 degrees of freedom. We assume that the upper parts will not move during motion. Thus we compute the joint trajectories only for the lower limbs. We have to plan the trajectories for the 12 joints of the legs. 3) Balance: The balance of humanoid robots can be defined thanks to the Zero Moment Point (ZMP). The ZMP is defined in [START_REF] Vukobratovic | Contribution to the synthesis of biped gait[END_REF] as a point, on the contact surface, where the total inertia force is equal to 0. If this point stays within the base of support, the robot keeps its balance. The position of the ZMP depends on the ground reaction effort.

ZMP s (t) ZMP f (t) = f (F(t)). (3)
ZMP s (t) and ZMP f (t) are the time history of the ZMP projected in the sagital and frontal planes and rely on the ground reaction effort, and hence on the joint trajectories {q(t), q(t), q(t)}.

B. Motion planning as SIP

The motion planning problem can be defined as a Semi-Infinite Programming (SIP) problem [START_REF] Vaz | Robot trajectory planning with semi-infinite programming[END_REF]. A SIP problem is an optimization problem with a finite number of variables and an infinite number of constraints [START_REF] Hettich | Semi-infinite programming: theory, methods, and applications[END_REF]. It consists in finding the parameter vector X that:

minimizes J( X,t) subject to ∀ j, ∀t ∈ [0, T ] g j ( X,t) ≤ 0 and ∀k h k ( X) = 0 (4)
Where F denotes the cost (or objective) function, g i the set of inequality constraint functions, h j the set of equality constraint functions.

1) Cost function: The choice of the cost function J(X,t) for motion planning must take into account the features of the robot and the desired application. For robot manipulators, some authors minimize motion duration [START_REF] Piazzi | Global minimum-time trajectory planning of mechanical manipulators using interval analysis[END_REF], or jerk [START_REF] Piazzi | Global minimum-jerk trajectory planning of robot manipulators[END_REF]. For humanoid robots, the energy consumption, taking into account the parameters of the motors (friction, ...) [START_REF] Miossec | Development of a software for motion optimization of robots -application to the kick motion of the hrp-2 robot[END_REF], or biological inspired, such as the minimum torque change [START_REF] Uno | Formation and control of optimal trajectory in human multijoint arm movement[END_REF], cost function can be considered. To improve the autonomy of the robots, we choose the cost function as the sum of the torque square:

F = n ∑ i=1 Γ 2 (t).dt (5)
2) Equality constraint functions: The set of the equality constraint functions h j (X) allows to define the motion. These functions are usually used to impose to constraints on some system state variables at given time instants such as the beginning or the end of a motion. For humanoid robot, we consider equality constraints as the position of the flying foot at the beginning, at the end and at mid-duration of the motion:

∀t k ∈ {0; T 2 ; T } h(X,t k ) = 0 (6)
3) Inequality constraint functions: The set of the inequality constraints g i (X) translates the physical limits of the system. For kicking motion we consider limits about joint position, velocity and torque and about sagital and frontal balance, as presented in Equation [START_REF]Safe motion planning computation for databasing balanced movement of humanoid robots[END_REF].

∀t ∈ [0, T ] ∀k ∈ {1, 2,...,n}

q k ≤ q k (t) ≤ q k qk ≤ qk (t) ≤ qk Γ k ≤ Γ k (t) ≤ Γ k ∀t ∈ [0, T ] ZMP s ≤ ZMP s (t) ≤ ZMP s ZMP f ≤ ŻMP f (t) ≤ ZMP f ( 7
)

JJ%

It is possible to take into account other constraint functions to avoid sliding or self-collision [START_REF] Escande | Continuous gradient proximity distance for humanoids free-collision optimized-postures[END_REF], for instance.

C. Time-Interval Discretization

In [START_REF]Safe motion planning computation for databasing balanced movement of humanoid robots[END_REF], we presented a new way for dealing with the inequality constraints : the time-interval discretization, which ensures constraint validity over whole motion duration and allows to use state-of-the-art algorithm such as IPOPT [START_REF]Introduction to IPOPT : a tutorial for downloading, installing and using IPOPT[END_REF]. This method uses Interval Analysis to ensure the inequality constraints validity over whole motion duration, by computing minimum and maximum values for the set of functions g i (t) when t is defined over a given interval [t].

Therefore the upper bound of g i (t): max g i are obtained in an easy and practical way by computing the upper bound of the inclusion function [g i ] for a time interval [t] [START_REF] Moore | Methods and Applications of Interval Analysis[END_REF]. The inequality constraints in (4) are replaced by:

∀i, ∀[t] ∈ IT Sup[g] i (X, [t]) ≤ 0 with IT = {[0,t 1 ], [t 1 ,t 2 ] ... [t k , T ]} (8)
In practice, the bounds thus derived may be too large because of over approximation in interval computing process. Still, there are several techniques that can be used to obtain tighter enclosures by using for instance Taylor series expansion or some global optimization techniques [START_REF] Hansen | Global optimization using interval analysis[END_REF].

D. Optimal kicking motion

We planned an optimal motion to make the robot kick a ball at the mid-duration of its motion at 0.6m/s speed. We impose that the collision will be at the position x = 1cm, h = 3cm as shown in Figure 3.

II. FAST MOTION RE-PLANNING

A. The Problem Under Study

In the previous section, we generated a kicking motion while assuming the location of the ball at x = 1cm. What happens if the ball is not at the expected position ? Figure 7 shows the result obtained when this optimal motion is used with the ball at x = 3cm. The foot hits the ball higher than expected. Thus the energy transmitted to the ball may be insufficient to reach the desired goal.

To improve the kicking motion, one solution could be to generate a new optimal motion with this new equality constraint. However, this solution is too time-consuming (about two hours for the previous one). In this sequel, we introduce a method which modifies the previous optimal kicking motion, in a very small CPU time while ensuring constraint satisfaction.

Our idea consists in replacing the set of inequality constraint ∀t ∈ [0, T ], g(X,t) < 0 by a set of bounds on the

B. Computation of the Feasible Sub-set

To make the robot able to adapt its motion to as many situations as possible we have to compute a feasible sub-set [X] that contains the optimal vector X and satisfy all the inequality constraint functions. Recent studies addressed the computation of feasible sets, using Interval Analysis, for the design of parallel or serial robots [START_REF] Ramdani | First results on the design of high speed parallel robots in presence of uncertainty[END_REF], [START_REF] Oetomo | Design strategy of serial manipulators with certified constraint satisfaction[END_REF]. In fact, we do not compute whole feasible set, but only an inner approximation of it. The sub-set [X] will be contained in the feasible set. We define a box as large as possible, then we solve the following problem:

maximize δ ∈ R + such as ∀i [X i ] = Xi + δ × [W i ] with ∀i 0 ∈ [W i ] ∀ j, ∀X ∈ [X], ∀t ∈ [0, T ] g j (X,t) < 0 ( 10
)
Where δ is the normalized width of the box and [W] a weighted interval vector that allows to ignore or give priority to some components of the box [X]. In this case, we propose that [W] is computing by using the distance between the optimal vector X and the first constraint violation along each direction.

As presented in Section I-C, we propose to use the timeinterval discretization which ensures constraint validity over whole motion duration. The inequality constraint in Equation ( 10) is replaced by:

∀i, ∀[t] ∈ IT Sup[g] i (X, [t]) ≤ 0 with IT = {[0,t 1 ], [t 1 ,t 2 ] ... [t k , T ]} (11)

C. Algorithm

The principle of the algorithm is to start from a large value of δ , and to reduce it by rejecting all the solutions in the corresponding box which violate a constraint.

Figure 8 shows the principle of this algorithm computing the feasible sub-set [X]. Using ALIAS software [START_REF]ALIAS-C++ A C++ Algorithms Library of Interval Analysis for equation Systems[END_REF], a branching algorithm with consistence tests, we search a box [z] that satisfies: with

[X] = X + δ k [W] find [z] ⊂ [X] such as ∃ j, ∃t ∈ [t] Sup[g] j ([z],t) > 0 (12)
Once the software finds a solution, [z], it stops and δ is chosen such that:

[z] ∩ X + δ k+1 [W] = / 0 ( 13
)
The algorithm stops when there is no solution to the problem [START_REF] Wieber | The humans toolbox, a homogeneous framework for motion capture, analysis and simulation[END_REF]. Eventually, computed X + δ f inal [W] is the largest box contained in the feasible set.

III. KICKING-MOTION ADAPTATION

A. Choice of the parameter to adapt

Obviously, it is not necessary to adapt all the motion parameters. Since, we are interesting in the collision location along the x-axis, we propose to adapt the trajectories of the knee, hip pitch and ankle pitch ( named LEG JOINT [START_REF] Piazzi | Global minimum-jerk trajectory planning of robot manipulators[END_REF][START_REF] Suleiman | On humanoid motion optimization[END_REF][START_REF] Lengagne | Guaranteed computation of constraints for safe path planning[END_REF] in Figure 2) which influence the motion in the sagital plane. The collision occurs at the half of the motion, thus we will only change the third B-splines parameters.

B. Feasible Sub-Set

Table I presents the result of our computation. Where the optimal value is Xi , the lower bound of the feasible set X i and the upper bound of the feasible set X i . The width of the feasible sub-set depends on the parameter, but it is interesting to see that some parameters can be changed within an interval of 5 degrees without making the robot falls ( since it ensures no constraint violation).

Figures 9 shows the time history of a joint value and the ZMP in the sagital plane for the optimal motion and their limits for all the motions contained in the feasible sub-set. In Figure 9 all the motions of the feasible sub-set will make the robot keeping its balance, since their limits stay within the size of the foot (-0.04 < ZMP s (t) < 0.068). Thus, we can choose any motion of the feasible sub-set without computing the inequality constraint functions to make sure that there is no constraint violation. All the motions of the feasible sub-set are safe. Figure 11 shows the feasible set of the couple (x, h) for all the motions in the feasible sub-set [X]. Unfortunately, it appears that we cannot achieve a kicking motion for the collision location (x = 3cm, h = 3cm). If we want the robot kicks the ball at 3cm high, the ball must be located between -1cm and 1.7cm. We choose to re-plan the optimal motion to make the robot kicks a ball at 3cm high and located at the position x = 1.5cm. Thus, we proceed to the optimization of the problem presented in Equation( 9) with these equality constraints:

C. Re-planned motion

find X ∈ [X] such as h( T 2 ) = 3cm x( T 2 ) = 1.5cm (14)
The optimization software spent less than one second of CPU time to find a solution. This solution is presented in Table II and the re-planned motion in Figure 12.

CONCLUSION

In this paper we presented the planning and fast replanning of safe motions. We applied our method to a kicking motion for a humanoid robot. The safe motion planning consists in solving a Semi Infinite Programming problem, using a time-interval discretization. Unfortunately, this method requires a large CPU time (2 hours), and hence cannot be done online. We generate a safe kicking motion which makes the HOAP-3 Humanoid Robot kick a ball located at 1cm from its foot. We showed that this motion is not good enough, if the ball is farther than expected. As a result we propose a safe re-planning method, which starting from the optimal motion, computes off-line a feasible subset of the motion parameters. By the way, we can achieve a fast re-planning which consists in finding, in this feasible sub-set, a solution that will validate a new set of equality constraint. In the future we will test this method for other motions, for example to adapt optimal step motions to a new position or to slopes.

Fig. 1 .

 1 Fig. 1. Time History of the Bsplines functions

Fig. 2 .

 2 Fig. 2. Denavit-Hartenberg parameters definition coordinate of the HOAP-3 Robot

Fig. 3 .

 3 Fig. 3. Scheme of the collision between the ball and foot of the robot

Fig. 4 .Fig. 5 .

 45 Fig. 4. Optimal kicking motion with a ball far from 1cm

Fig. 6 .Fig. 7 .

 67 Fig. 6. Time history of the ZMP in the frontal plane

Fig. 8 .

 8 Fig. 8. Example of a feasible set and of its inner approximation : the feasible sub-set [X]

Fig. 9 .

 9 Fig. 12. Re-planned kicking motion

JJA

Fig. 10 .

 10 Fig.10. Time history for the x-position of the flying foot for the optimal motion, for the re-planned motion, and its limits for all the solution in the feasible sub-set

TABLE I TABLE OF THE

 IOF (b), the ZMP in the sagital plane is presented, VALUE FOR THE OPTIMAL MOTION, FOR THE FEASIBLE SUB-SET (ANGLE ARE GIVEN IN DEGREE)

	position	Xi	X i	X i
	Right Hip Pitch	-22.60	-26.58	-20.34
	Right Knee	39.68	38.20	41.21
	Right Ankle Pitch	-17.33	-18.88	-16.05
	Left Hip Pitch	-17.12	-19.13	-13.97
	Left Knee	19.56	13.97	22.90
	Left Ankle Pitch	4.51	-3.66	5.90
				JJJ

TABLE II TABLE OF

 IIOF

THE VALUE FOR THE OPTIMAL MOTION AND FOR THE RE-PLANNED MOTION (ANGLE ARE GIVEN IN DEGREE)