Left-handed optical radiation torque

Davit Hakobyan and Etienne Brasselet * Optical forces and torques are two mechanical degrees of freedom available to manipulate matter, and form the basis of optical tweezing strategies [START_REF] Grier | A revolution in optical manipulation[END_REF][START_REF] Padgett | Tweezers with a twist[END_REF] . In contrast to the Keplerian intuition that objects should be pushed downstream an incident photon flux, the concept of 'negative' optical forces has recently been described [START_REF] Lee | Optical solenoid beams[END_REF][START_REF] Sukhov | On the concept of 'tractor beams[END_REF] and has triggered many developments [START_REF] Salandrino | Reverse optical forces in negative index dielectric waveguide arrays[END_REF][START_REF] Chen | Optical pulling force[END_REF][START_REF] Sukhov | Negative nonconservative forces: optical "tractor beams" for arbitrary objects[END_REF][START_REF] Novitsky | Single gradientless light beam drags particles as tractor beams[END_REF][START_REF] Nemirovsky | Negative radiation pressure and negative effective refractive index via dielectric birefringence[END_REF][START_REF] Novitsky | Material-independent and size-independent tractor beams for dipole objects[END_REF][START_REF] Dogariu | Optically induced 'negative forces[END_REF][START_REF] Brzobohatý | Experimental demonstration of optical transport, sorting and self-arrangement using a 'tractor beam[END_REF][START_REF] Wang | Dynamical and phase-diagram study on stable optical pulling force in Bessel beams[END_REF][START_REF] Kajorndejnukul | Linear momentum increase and negative optical forces at dielectric interface[END_REF] . Here, we report on the counterintuitive angular analogue of negative optical forces by demonstrating that circularly polarized Gaussian light beams give rise to torque with opposite sign to that of the incident optical angular momentum. Such a 'left-handed' mechanical effect is demonstrated by the use of an inhomogeneous and anisotropic transparent macroscopic medium. Practical difficulties associated with the direct observation of optically induced spinning of a macroscopic object are circumvented via the rotational Doppler effect [START_REF] Allen | A radiation torque experiment[END_REF][START_REF] Garetz | Angular Doppler effect[END_REF] . These results shed light on spin-orbit optomechanics and equip the left-handed optomechanical toolbox with angular features.

The redirection of a photon flux by a material system is essen tially associated with non conservative optical radiation forces. In particular, negative optical forces occur when the linear momentum of light experiences a net forward scattering that depends on the structure of the incident electromagnetic fields, the properties of the objects and their surrounding environment [START_REF] Dogariu | Optically induced 'negative forces[END_REF] . Quite naturally, because scattering is a universal feature common to all waves, the concept of negative forces is also found in acoustics [START_REF] Marston | Axial radiation force of a Bessel beam on a sphere and direction reversal of the force[END_REF][START_REF] Marston | Negative axial radiation forces on solid spheres and shells in a Bessel beam[END_REF][START_REF] Mitri | Negative axial radiation force on a fluid and elastic spheres illuminated by a high-order Bessel beam of progressive waves[END_REF] . Intriguingly, its angular counterpart has not yet emerged. This corresponds to a light field exerting an optical radiation torque (ORT) on matter, where the direction of the torque is reversed with respect to that of the incident angular momentum assumed to be at play to produce the ORT. We shall refer to such a reverse angular optomechanical effect as 'left handed' instead of 'negative' to avoid confusion, as it may occur either with τ z > 0 or τ z < 0, where τ z is the projection of the ORT along z that defines the pro pagation direction of the light. Accordingly, we will refer to 'right handed' for the intuitive common situation. Here, we report the observation of left handed ORT for the case of a circularly polarized Gaussian light beam carrying an angular momentum per photon of σħ along its propagation direction z (σ = ±1, ħ is the reduced Planck constant), namely στ z < 0, whatever σ.

The non trivial nature of left handed ORT is exemplified by recalling the two types of light matter interaction processes leading to ORT that are driven by spin angular momentum. The first is the dissipative case where the exchange of angular momen tum is driven by absorption, as experimentally reported in ref. 20. By considering an absorbing optically isotropic object characterized by an overall power attenuation coefficient of 0 ≤ ε ≤ 1, the angular momentum balance gives τ z = σεħ (Fig. 1a). Because στ z = εħ ≥ 0, left handed ORT cannot occur. The second case refers to non dissipative exchange of angular momentum involving a birefringent medium, as demonstrated in 1936 using a centimetre sized crystalline half wave plate [START_REF] Beth | Mechanical detection and measurement of the angular momentum of light[END_REF] , inspiring significant advances in down sizing down to the micrometre scale [START_REF] Nieminen | Optical alignment and spinning of laser-trapped microscopic particles[END_REF] . The angular momentum balance gives τ z = σ(1 -cosΔ)ħ, where Δ is the birefringent phase retardation of the wave plate (Fig. 1b). As στ z = (1 -cosΔ)ħ ≥ 0, left handed ORT is also prohibited.

Our idea to achieve left handed ORT consists of coupling the incident spin angular momentum to the spatial degrees of freedom of the light field, namely orbital angular momentum. Indeed, if ORT is dominated by the orbital contribution rather than that of the spin, reversing the angular momentum balance becomes possible. Because such spin orbit light scattering basically takes place in optically inhomogeneous and anisotropic transparent media, we propose, for the purpose of demonstration, to restrict the present study to space variant uniaxial slabs with constant birefrin gent phase retardation Δ and optical axis orientation angle ψ of the form ψ(ϕ) = qϕ, where q is the half integer and ϕ the polar angle in the plane of the slab. When Δ = π, a collimated circularly polarized Gaussian beam carrying an incident spin angular momentum per photon s in z = σhthat impinges on axis onto the azimuthally pat terned slab will emerge from it with a total angular momentum per photon of j out z = s out z + l out z , where s out z = -σhand l out z = 2σqhare respectively the spin and orbital contributions [START_REF] Marrucci | Optical spin-to-orbital angular momentum conversion in inhomogeneous anisotropic media[END_REF] . The angular momentum balance τ z = j in zj out z thus gives

τ z = 2σ 1 -q h - (1)
Because στ z = 2(1q)ħ whatever σ, left handed ORT is expected when q ≥ 3/2, whereas right handed ORT takes place when q ≤ 1/2, as shown in Fig. 1c. In addition, the particular case of q = 1 corresponds to vanishing ORT, as expected from the rotational invariance around z.

In practice we used azimuthally patterned birefringent disk shaped glass slabs with q = (1/2, 1, 3/2) and Δ = π at a wavelength of 532 nm. The spatial distribution of the slow axis is depicted in Fig. 2 (first column). The samples were fabricated by femtosecond laser writing of self assembled nanogratings in silica glass [START_REF] Shimotsuma | Self-organized nanogratings in glass irradiated by ultrashort light pulses[END_REF] , and the induced form birefringence patterns were achieved by control ling the writing parameters, in particular the polarization azimuth of the writing beam [START_REF] Beresna | Radially polarized optical vortex converter created by femtosecond laser nanostructuring of glass[END_REF] . Their polariscopic optical characterization at a wavelength of 532 nm is summarized in Fig. 2, where crossed linear polarizer imaging (second column) emphasizes the azimuthal patterning of the optical axis, and crossed and parallel circular polarizer imaging (third and fourth columns) present the uniform birefringent phase retardation Δ = π. See also the 'Materials and Methods' section.

A straightforward experimental approach for observing the optical torque would consist of direct observation of the light induced spinning of a free to rotate sample around the z axis, the left handed nature of the ORT being retrieved from a sense of rotation opposite that prescribed by the incident spin angular momentum. However, considering a transfer of ħ angular momen tum per photon, a light beam with power P and angular frequency ω produces a torque of magnitude Γ = P/ω, so Γ ≈ 1 × 10 -16 N m University of Bordeaux, CNRS, Laboratoire Ondes et Matière d'Aquitaine, 351 cours de la libération, F-33400 Talence, France. *e-mail: e.brasselet@loma.u-bordeaux1.fr under 1 W visible light illumination. Then, assuming free rotation in an inviscid medium one finds a rotation frequency v that depends on the illumination time T following v(T ) = ΓT/(πmR 2 ), where m is the mass of the slab. According to the sample parameters (the sample is structured over an area of diameter d = 4 mm centred on a disk shaped glass substrate of radius R = 12.7 mm, thickness H = 3 mm and density ρ = 2.5 × 10 3 kg m -3 ) one typically obtains v = 1 mHz after a three month illumination at P = 1 W, which clearly prevents a realistic experimental implementation.

The practical difficulty of observing the mechanical consequence of the spin orbit scattering process on matter is circumvented by probing the mechanical effect of the ORT on the light itself. This is implemented using a rotational Doppler experiment [START_REF] Allen | A radiation torque experiment[END_REF][START_REF] Garetz | Angular Doppler effect[END_REF] by rotat ing the slab at controlled angular velocity Ω around the z axis. Non zero ORT therefore causes a frequency shift δω for the output light, whose expression can be derived from energy conservation, namely w + ħδω = 0 where w = τ z Ω is the work per photon produced on the slab. This gives

δω = 2σΩ q -1 (2)
Because τ z = ħδω/Ω, we stress that the determination of the rotational Doppler frequency shift is actually an ORT measurement.

The experiment was carried out by observing the intensity pat terns that result from the non collinear superposition of the output light field and a reference Gaussian beam with angular fre quencies ω + δω and ω, respectively, as outlined in Fig. 3a. Indeed, we observe 2q fork interference patterns that are reminiscent of the optical phase singularity with topological charge 2qσ emerging from the sample [START_REF] Marrucci | Optical spin-to-orbital angular momentum conversion in inhomogeneous anisotropic media[END_REF] , as illustrated by the snapshots shown in Fig. 3b d for q = (1/2, 1, 3/2). The determination of δω is achieved by analysing the dynamics of these intensity patterns, denoted I(x, y, t), with time t. Indeed, by constructing spatiotem poral interferograms I(x = x 0 , y, t) we obtain fringing patterns that may or may not drift along the y axis depending on q, σ and Ω. This is illustrated in Fig. 3e g, for q = (1/2, 1, 3/2) when σΩ > 0.

q = 1 q = 3/2 q = 1/2 P in P out
The sign of the latter drift gives access to the sign of δω from a mere visual inspection. In other words, with our set up, when σΩ > 0, a drift towards y < 0 refers to δω < 0 (Fig. 3e), whereas a drift towards y > 0 indicates δω > 0 (Fig. 3g), and no drift implies δω = 0 (Fig. 3f) (Supplementary Section 1). Because στ z = ħδω/σΩ, our results indicate that στ z < 0 when q = 3/2, hence demonstrating that left handed ORT takes place.

Regarding the magnitude of δω, we define the time dependent quantity I(t) = ∫ I(x 0 , y, t)dy -〈∫ I(x 0 , y, t)dy〉 t (Fig. 4a), where 〈•〉 t holds for time averaging, where the spectrum exhibits a well defined peak at frequency |δω|/2π (Fig. 4b), as demonstrated in Supplementary Section 2. The dependence of δω as a function of Ω is shown in Fig. 4c,d, corresponding to the range |δω|/ω = 10 -17 10 -16 for q = (1/2, 3/2) and σ = ±1. In Fig. 4c,d, the solid lines refer to the theoretical results given by equation ( 2), from which we can infer the excellent agreement between the experimen tal data and theoretical predictions. We thus obtain the precise experimental determination of both the sign and magnitude of the rotational Doppler frequency shift, thereby allowing the quanti tative identification of left handed ORT.

Next we address the question of whether left handed ORT may occur for Δ ≠ π. In this case the incident light field emerges from the slab with a total angular momentum per photon of j out z = σ[cosΔ + 2qsin 2 (Δ / 2)]h -(Supplementary Section 3) from which results an ORT of

τ z = 2σ 1 -q sin 2 Δ/2 ( )h - (3)
Previous conclusions are therefore unaltered whatever the value of Δ, up to a factor of sin 2 (Δ/2) for the ORT magnitude.

Experimentally, this situation is explored by detuning the incident wavelength from 532 nm, setting it at 632.8 nm, so Δ = (532/632.8)π (that is, 0.84π). The results are displayed in Fig. 5, where the rotational Doppler effect is assessed independently for the two output light field components with opposite helicities ±σ and q = (1/2, 1, 3/2). This is done by using a circularly polarized reference beam with helicity σ ref = σ (Fig. 5b d) and σ ref = σ (Fig. 5e g), respectively. We conclude from Fig. 5b d that photons whose helicity is unchanged did not produce ORT. In con trast, helicity flipped photons give rise to right handed, zero or left handed ORT depending on q, as shown in Fig. 5e g. Noticeably, the factor f = sin 2 (Δ/2) in equation ( 3) corresponds to the fraction of output photons that produce work on the sample (Supplementary Section 3). In addition, because the energy conservation condition per photon should be written 'per working photon', hence τ z Ω + fħδω = 0, the predicted magnitude of δω is unchanged with respect to the case Δ = π. This has been quantitatively verified and can be qualitatively grasped from a comparison of Fig. 3e g and Fig. 5e g. We thus generalize the proposed concept of left handed ORT to arbitrary values of birefringent phase retardation. The above results remain valid for any superposition of left and right handed circularly polarized Gaussian beams. Indeed, account ing for an electric field of the form αG +1 + βG -1 where (α, β) ∈ C 2 and G σ refers to a Gaussian field with helicity σ, one obtains

s in z τ z = 2[(|α| 2 -|β| 2 ) / (|α| 2 + |β| 2 )] 2 (1 -q)sin 2 (Δ/2) h -2
, which leaves the previous conclusions unchanged. Namely, s in z τ z < 0 for q > 1 whatever s in z .

In the context of previous experimental studies on light induced rotation of optically isotropic transparent non axisymmetric micro objects, we note that a Gaussian beam has been shown to lead to either clockwise only or anticlockwise only spinning motion inde pendently of the incident angular momentum [START_REF] Higurashi | Optically induced rotation of anisotropic micro-objects fabricated by surface micromachining[END_REF][START_REF] Galajda | Rotation of microscopic propellers in laser tweezers[END_REF] , and therefore cannot be considered as a left handed optomechanical manifestation. This is also the case for the reversed orbiting motion of a micropar ticle with irregular shape induced by a linearly polarized vortex beam (reported in ref. 28), because the observed sense of rotation is inde pendent of the sign of the incident angular momentum, and also for the polarization dependent rotational motion of an achiral micro particle with well controlled shape induced by a polarized vortex beam (reported in ref. 29), where the observed sense of rotation is given by the sign of the incident angular momentum, hence corre sponding to a right handed ORT. On the other hand, in the case of optically induced torque driven by absorption only, the sign of the The sample is fixed on a rotating stage whose angular velocity Ω is externally controlled. The output light field, labelled E out (ω', σ), has frequency ω' = ω + δω, which may differ from ω, and helicity opposite to that of the input. b-d, Snapshots of the 2q-fork intensity patterns I(x, y, t i ) at time t = t i between the output beam and a noncollinear reference beam with helicity σ ref = σ and frequency ω for q = 1/2, q = 1 and q = 3/2, respectively. e-g, Spatiotemporal interferograms taken at x = x 0 , I(x 0 , y, t), when σΩ > 0 for q = 1/2 (e), q = 1 (f) and q = 3/2 (g). The drift direction of the fringing patterns allows the determination of the sign of the rotational Doppler frequency shift, namely δω < 0, δω = 0 or δω > 0. The wavelength is 532 nm and Ω = 0.37 rad s 1 .

E re f ( ω ,σ re f ) E in ( ω ,σ ) E o u t ( ω ', -σ ) Δ = π Ω q = 1/2 b c d e f g q = 1 q = 3/2 x = x 0 t = t i
torque is always the same as that of the incident angular momentum, whatever the polarization state and orbital content of the incident light field (shown in ref. 30), which precludes a left handed effect. In the general case where absorption and/or birefringence and/or scattering can be at play, the circumstances under which a left handed rotational effect takes place certainly deserve further study. By addressing the rotational degree of freedom, the present study strengthens the emergence of a novel field in optical trapping and micromanipulation left handed optomechanics driven by spin orbit scattering of light. As already emphasized for the case of negative optical forces, a left handed optomechanical toolbox may find application in the metrology of fluidic and soft matter systems [START_REF] Dogariu | Optically induced 'negative forces[END_REF][START_REF] Brzobohatý | Experimental demonstration of optical transport, sorting and self-arrangement using a 'tractor beam[END_REF][START_REF] Kajorndejnukul | Linear momentum increase and negative optical forces at dielectric interface[END_REF] . Original optical micro/nanostructuring of matter could also be envisioned for photonics applications. More generally, this study addresses the interaction of light with structured media and is likely to lead to new developments in the science and appli cations of structured light.

Methods

The three samples with q = (1/2, 1, 3/2) were form-birefringent nanostructured glass slabs from Altechna R&D. The azimuthal distribution of their optical axis is shown to satisfy ψ = qϕ, where ϕ is the polar angle in the plane of the slab. The glass substrates had a diameter of 1 inch, a thickness of 3 mm, and the structured area corresponded to a 4-mm-diameter region centred on the substrate. The samples were designed to exhibit a uniform birefringent phase retardation over the whole structured area, which equalled π at a wavelength of 532 nm. The polariscopic characterization of the samples was carried out using incoherent illumination at 532 nm provided by a halogen lamp, in front of which were placed tracing paper and an interference filter (wavelength of 532 nm) with 10 nm spectral width (Thorlabs). The rotational Doppler experiments were performed using either a diode-pumped solid-state laser (Cobolt) operating at 532 nm or a He Ne gas laser (Melles Griot) operating at 632.8 nm. Both laser sources operated on the fundamental TEM 00 Gaussian mode. The incident light beam impinged at normal incidence onto the slab (and was centred on it), and was collimated into a ∼2-mm-diameter beam (at 1/e 2 of its maximum intensity) using a homemade optical telescope made of a pair of planoconvex lenses whose focal lenses were adapted to the laser source used. The slabs were placed in a rotating holder whose angular velocity was controlled by a voltagecontrolled d.c. motor. The interference intensity patterns were recorded using a CMOS video camera (Thorlabs) and the videos were processed using homemade Matlab code. In all experiments the polarization state of the light was ensured using either linear or circular polarizers (Meadowlark Optics). , q = 1 (c) and q = 3/2 (d), and in the case σσ ref < 0 for q = 1/2 (e), q = 1 (f) and q = 3/2 (g). Here, the wavelength is 632.8 nm, so Δ = 0.84π and Ω = 0.37 rad s 1 .

Figure 2 |Figure 1 |

 21 Figure 2 | Spin-orbit light scattering samples. a-c, Slow-axis spatial distribution of the form-birefringent nanostructured glass plates with q = 1/2, 1 and 3/2, where segments indicate local orientation. d-l, Polariscopic analysis performed by optical imaging between crossed linear polarizers (d-f), crossed circular polarizers (g-i) and parallel circular polarizers (j-l) under incoherent spatially extended illumination at a wavelength of 532 nm. P in and P out denote the input and output polarization states of light. The structured area has a circular cross-section with diameter d = 4 mm.

Figure 3 |

 3 Figure 3 | Demonstration of right-handed, zero and left-handed ORT. a, Schematic of the experiment. A circularly polarized light beam with helicity σ = ±1 and angular frequency ω, labelled E in (ω,σ), is normally incident on a sample designed to have a uniform birefringent phase retardation Δ = π at a wavelength of 532 nm.The sample is fixed on a rotating stage whose angular velocity Ω is externally controlled. The output light field, labelled E out (ω', σ), has frequency ω' = ω + δω, which may differ from ω, and helicity opposite to that of the input. b-d, Snapshots of the 2q-fork intensity patterns I(x, y, t i ) at time t = t i between the output beam and a noncollinear reference beam with helicity σ ref = σ and frequency ω for q = 1/2, q = 1 and q = 3/2, respectively. e-g, Spatiotemporal interferograms taken at x = x 0 , I(x 0 , y, t), when σΩ > 0 for q = 1/2 (e), q = 1 (f) and q = 3/2 (g). The drift direction of the fringing patterns allows the determination of the sign of the rotational Doppler frequency shift, namely δω < 0, δω = 0 or δω > 0. The wavelength is 532 nm and Ω = 0.37 rad s 1 .

δω (rad s - 1 )Figure 4 |

 14 Figure 4 | Quantitative analysis of rotational Doppler experiments. a, Time-oscillating signal I = ∫ I(x 0 , y, t)dy -〈∫ I(x 0 , y, t)dy〉 t , where I refers to the intensity interferogram pattern shown in Fig.3and 〈•〉 t refers to time averaging. b, Fast Fourier transform (FFT) spectrum of Ī, which exhibits a well-defined peak at frequency |δω|/2π. Here Ω = 0.37 rad s 1 and the wavelength is 532 nm. c,d, Rotational Doppler angular frequency shift δω = ω' ω as a function of the rotating slab angular velocity Ω for q = 1/2 and q = 3/2 when σ = 1 (c) and σ = +1 (d). Solid lines show the theoretical results given by equation (2).

Figure 5 |

 5 Figure5| Generalized ORT rotational Doppler experiments. a, Schematic of the experiment when the uniform birefringent phase delay of the nanostructured glass slab is Δ(ω) ≠ π. This leads to an elliptically polarized output light beam that consists of the coherent superposition of two fields with opposite helicities, E out (ω, σ) and E out (ω', σ), whose angular frequencies and relative intensity weights may differ. b-g, The contributions of the output light field components to ORT are retrieved from their respective spatiotemporal interferograms using a non-collinear reference beam with helicity σ ref = ±σ and frequency ω, here shown for σΩ > 0: in the case σσ ref > 0 for q = 1/2 (b), q = 1 (c) and q = 3/2 (d), and in the case σσ ref < 0 for q = 1/2 (e), q = 1 (f) and q = 3/2 (g). Here, the wavelength is 632.8 nm, so Δ = 0.84π and Ω = 0.37 rad s 1 .

Acknowledgements

The authors thank C. Loussert for the polariscopic optical characterization of the samples. This study received financial support from the French State in the frame of the 'Investments for the future' Programme IdEx Bordeaux (reference ANR-10-IDEX-03-02).

Author contributions


S Ébastien

Aline Kurtzmann

SOME PARTICULAR SELF-INTERACTING DIFFUSIONS: ERGODIC BEHAVIOUR AND ALMOST SURE CONVERGENCE

Keywords: 60K35, 37C50 Self-interaction diffusion, Reinforced processes, Stochastic approximation

This paper deals with some self-interacting diffusions (Xt, t ≥ 0) living on R d . These diffusions are solutions to stochastic differential equations: dXt = dBt -g(t)∇V (Xt -µ t )dt, where µ t is the empirical mean of the process X, V is an asymptotically strictly convex potential and g is a given function. We study the ergodic behaviour of X and prove that it is strongly related to g. Actually, we show that X is ergodic (in the limit-quotient sense) if and only if µ t converges a.s. We also give some conditions (on g and V ) for the almost sure convergence of X.

Introduction

Processes with path-interaction have been an intensive research area since the seminal work of Norris, Rogers and Williams [START_REF] Norris | Self-avoiding random walk: a Brownian motion model with local time drift[END_REF]. More precisely, selfinteracting diffusions have been first introduced by Durrett and Rogers [START_REF] Durrett | Asymptotic behaviour of Brownian polymers[END_REF] under the name of Brownian polymers. They proposed a model for the shape of a growing polymer. Denoting by X t the location of the end of the polymer at time t, X satisfies a SDE with a drift term depending on its own occupation measure (in dimension 1, we define it through the local time of X). One is then interested in rescaling X (see [START_REF] Cranston | Self-attracting diffusions: two cases studies[END_REF][START_REF] Cranston | The strong law of large numbers for a Brownian polymer[END_REF][START_REF] Roynette | Boundedness and convergence of some selfinteracting diffusions[END_REF][START_REF] Mountford | An asymptotic result for Brownian polymers[END_REF][START_REF] Raimond | Self-attracting diffusions: case of constant interaction[END_REF]). Later, another model of polymers has been proposed by Benaïm, Ledoux and Raimond [START_REF] Benaïm | Self-interacting diffusions[END_REF]. They have studied a class of self-interacting diffusions depending on the empirical measure. When the process is living on a compact Riemannian manifold, they have proved that the asymptotic behaviour of the empirical measure can be related to the analysis of some deterministic dynamical flow defined on the space of the Borel probability measures. Benaïm and Raimond [START_REF] Benaïm | Self-interacting diffusions III: symmetric interactions[END_REF] went further in this study and in particular, they gave sufficient conditions for the a.s. convergence of the empirical measure. Very recently, Raimond [START_REF] Raimond | Self-interacting diffusions: a simulated annealing version[END_REF] has generalized the previous work: he has studied the asymptotic properties of a process X, living on a Riemannian compact manifold M , solution to the SDE dX t = dB t -g(t)∇V * µ t (X t )dt, (

with V * µ t (x) = 1 t t 0 V (x, X s )ds, µ t = 1 t t 0 δ Xs ds, |g(t)| ≤ a log(t) and g (t) = O(t -γ ) with 0 < γ ≤ 1. He has proved that, unless g is constant, the approximation of µ t by a deterministic flow is no longer valid. He has more particularly investigated the example M = S n and V (x, y) = -cos d(x, y) (where d is the geodesic distance on S n ) and proved that a.s. µ t converges weakly towards a Dirac measure. For an overview on reinforced processes, we refer the reader to Pemantle's survey [START_REF] Pemantle | A survey of random processes with reinforcement[END_REF].

In the present paper, we are concerned with some self-interacting processes living on R d . Consider a smooth potential V : R d → R + and an application g : R + → R * + . Our goal is to study the ergodic behaviour of the self-interacting diffusion X solution to dX t = dB t -g(t)∇V (X t -µ t )dt, X 0 = x.

(1.2)

where B is a standard Brownian motion and µ t denotes the empirical mean of X:

µ t = 1 r + t r μ + t 0 X s ds , µ 0 = µ. (1.3)
Here µ is an initial (given) probability measure on R d , μ denotes the mean of µ and r > 0 is an initial weight (it permits to consider any initial probability measure). First, note that for a quadratic interaction potential V , the process satisfying (1.2) is exactly of the form of (1.1) and, in both cases, the occupation measure is penalized by g(t). Afterwards, a natural generalization of this process is the class of self-interacting diffusions discussed here. The interesting point is that we manage to study precisely the asymptotic behaviour of X and prove a convergence criterion. Moreover, this model could be used to represent the behaviour of social insects, as the ants trails. Indeed, ants mark their paths with the trails' pheromones. Certain ants lay down an initial trail of pheromones as they return to the nest with food. This trail attracts other ants and serves as a guide. As long as the food source remains, the pheromone trail will be continually renewed. Despite the quick evaporation, the path is reinforced and so, the ants manage to gradually find the best route. In this (simplified) model, the function g reflects the speed of evaporation and X denotes the trail.

In order to study the solution to (1.2), it is natural to introduce the process Y , defined by

Y t = X t -µ t . (1.4) It is easily seen that (Y t , t ≥ 0) is the solution to the SDE dY t = dB t -g(t)∇V (Y t )dt -Y t dt r + t , Y 0 = x -µ; (1.5) and dµ t = Y t dt r+t .
As Y is a (non-homogeneous) Markov process, it is easier to study Y than X. Indeed, we will prove that Y converges a.s. and satisfies the pointwise ergodic theorem. Due to that, the behaviour of X could seem a bit easy at first glance. But it really shows unexpected behaviours and, in particular, it does not satisfy the pointwise ergodic theorem in general (because µ t does not converge, except for functions g going fast to infinity). This explains how difficult is the study of more general self-interacting diffusions in non-compact spaces (see Kurtzmann [START_REF] Kurtzmann | The ODE method for some self-interacting diffusions on R d[END_REF]), driven by the generic equation dX t = dB t -R d ∇V (X t , x)dµ t (x)dt. The remainder of the paper is organized as follows. First, we enumerate the hypotheses and state the main results in Section 2. We motivate our study, in Section 3, by the basic case when V is quadratic, for which we have an explicit expression of X (in terms of Brownian martingale). Section 4 deals with the description of the behaviour of Y near the local extrema of V . Finally, Section 5 is devoted to the proof of the main results.

Technical assumptions and main results

In the sequel, (•, •) stands for the Euclidian scalar product. We also denote by P(R d ) the set of probability measures on R d .

Consider the potential V : R d → R + . Let M ax = {M 1 , . . . , M p } be the (finite) set of saddle points and local maxima of V and denote by M in = {m 1 , . . . , m n } the (finite) set of the local minima of V . So M in ∪ M ax is the set of critical points of V . We assume that V is either quadratic (Section 3) or:

1) (regularity and positivity) V ∈ C 2 (R d ) and V > 0;

2) (convexity) V = W + χ where χ is a compactly supported function such that ∇χ is C-Lipschitz (with C > 0) and there exists c > 0 such that ∇ 2 W ≥ cId; 3) (growth) there exists a > 0 such that for all x ∈ R d , we have ∆V (x) ≤ aV (x) and lim

|x|→∞ |∇V (x)| 2 V (x) = ∞; (2.1) 4) (critical points) ∀m i , ∀ξ ∈ R d , (∇ 2 V (m i )ξ, ξ) > 0 and for all M i , ∇ 2 V (M i ) admits a negative eigenvalue. Remark 2.1. By the growth condition (2.1), |∇V | 2 -∆V is bounded by below.
Suppose also that g : R + → R + is nondecreasing and g ∈ C 1 (R + ). We denote by g(∞) the limit of g(t) and we exclude the trivial case where g is identically zero, so that g(∞) > 0. Let G(t) := t 0 g(s)ds and G -1 be its generalized inverse:

G -1 (t) := inf{u ≥ 0; G(u) ≥ t}. Remark 2.2. If g(∞) = ∞, then for all T > 0, we have that G -1 (t + T ) - G -1 (t) -→ t→∞ 0.
The following easy result will be very useful in the sequel.

Lemma 2.3. Suppose that g (t)/g 2 (t) converges to 0. Then the following hold: i) for any c > 0, t 0 s 2 e 2cG(s) ds = O(t 2 e 2cG(t) /g(t)); ii) if g(∞) = ∞, then we have

t 0 g (s) g(s) 2 G(s)ds = O(t); iii) for H(t) := t 0 e -cG(u)
(r+u) 2 du, the following expansion holds:

H(t) = H(∞) - 1 cg(t)(r + t) 2 e -cG(t) + o e -cG(t) t 2 g(t)
.

Proof. We deduce all these estimates from an integration by parts:

t 0 s 2 e 2cG(s) ds = t 2 e 2cG(t) 2cg(t) - t 0 s g(s) - s 2 g (s) 2g(s) 2 e 2cG(s) c ds = O(t 2 e 2cG(t) /g(t)),
and we obtain H(t)-H(s) = e -cG(s) r+s -e -cG(t) r+t -c t s g(u)e -cG(u) du r+u . Similarly for t large enough and u such that g(u) > 0, we find

t u g (s) g(s) 2 G(s)ds = -G(t) g(t) + G(u) g(u) + t -u = O(t).
2.1. Existence. We begin by proving that the SDE studied admits a unique global strong solution.

Proposition 2.4. For any x ∈ R d , µ ∈ P(R d ) and r > 0, there exists a unique global strong solution (X t , t ≥ 0) of (1.2).

Proof. The local existence and uniqueness of the solution to (1.2) is standard. We only need to prove here that Y , hence X (since X t := Y t + t 0 Y s ds r+s ), does not explode in a finite time. To this aim, apply Itô's formula to the function x → V (x):

dV (Y t ) = (∇V (Y t ), dB t ) + 1 2 ∆V (Y t ) -g(t)|∇V (Y t )| 2 - 1 r + t (∇V (Y t ), Y t ) dt,
and introduce the sequence of stopping times τ 0 = 0 and

τ n = inf{t ≥ 0; V (Y t ) + t 0 g(s)|∇V (Y s )| 2 ds > n}.
By --→ f dµ ∞ . Theorem 2.6.

(1) The process Y satisfies the pointwise ergodic theorem: almost surely, the empirical measure of Y converges weakly to a measure, which is a convex combination of Dirac measures taken in the local minima of V .

(2) The process X satisfies the pointwise ergodic theorem if and only if the mean-process µ t converges almost surely.

A necessary condition for the convergence of µ t is that 0 is the unique minimum of V . We will prove this result in §5.

A motivating example

We consider V (x) = 1 2 (x, cx), where c is a symmetric positive definite matrix. Let X be the solution of the SDE

dX t = dB t -g(t) cX t - r r + t cµ - 1 r + t t 0 cX s ds dt, X 0 = x. (3.1)
Without any loss of generality, we suppose that d = 1. When d ≥ 1, it suffices to diagonalize the matrix c and to remark that, for an orthogonal matrix U , the process (U • B s , s ≥ 0) is also a Brownian motion.

Explicit expression of X.

Lemma 3.1. If X is the solution to (3.1), then we have

Y t := X t -μt = e -cG(t) r + t t 0 (r + s)e cG(s) dB s + r(x -µ) .
Proof. The process Y satisfies

dY t = dB t -cg(t) + 1 r + t Y t dt, Y 0 = x -µ. (3.2)
To express Y in terms of a Brownian martingale, we consider the function of Y defined by U t := (r + t)e cG(t) Y t . Then Itô's formula implies

dU t = (r + t)e cG(t) dB t , U 0 = r(x -µ). Corollary 3.2. Let F (t) = t 0 e -cG(s) g(s)
r+s ds. The solution to the SDE (3.1) is given by

X t = x + rc(µ -x)F (t) + t 0 1 -(r + s)ce cG(s) (F (t) -F (s)) dB s .
Proof. Remark that d dt µ t = Yt r+t . So, by Fubini's theorem for stochastic integrals (see [START_REF] Revuz | Continuous Martingales and Brownian Motion[END_REF] p.175), we have

µ t = t 0 (r + s)e cG(s) (H(t) -H(s))dB s + r(x -μ)H(t) + μ with H(t) = t 0 e -cG(u) (r+u) 2 du = 1 r -cF (t) -e -cG(t)
r+t . As X t = Y t + µ t , the latter result implies the desired expression.

3.2. Ergodic result. We begin to prove the pointwise ergodic theorem for the following non-homogeneous (Gauss-)Markov process.

Proposition 3.3. Let a : R + → R + be a continuous function, A(t) := t 0 a(s)ds and K(t) := e -2A(t) t 0 e 2A(s) ds. Suppose that a(∞) = lim t→∞ a(t)
exists and is nonzero, so that

K(∞) = 1 2a(∞) < ∞. Consider the process dZ t = -a(t)Z t dt + dB t , Z 0 = z.
Then, denoting by γ the centered Gaussian measure with variance K(∞) (with the convention γ = δ 0 for K(∞) = 0), we have for all continuous bounded function ϕ

1 t t 0 ϕ(Z s )ds a.s.
--→ t→∞ ϕ(z)γ(dz).

Proof. We prove the result for the Fourier transform. First, note that

Z t = e -A(t) t 0 e A(s) dB s + z . Let F s := σ(B u , 0 ≤ u ≤ s). Knowing F s , Z t has the Gaussian dis- tribution with mean m(s, t) := e -(A(t)-A(s)) Z s and variance K(s, t) := e -2A(t) t s e 2A(u) du. Fix t ∈ R + , u ∈ R and define the martingale M t,u s := E e iuZt |F s = exp ium(s, t) -u 2 2 K(s, t) . Applying Itô's formula to s → M t,u
s , we find that dM t,u s = iue -(A(t)-A(s)) M t,u s dB s , and so

e iuZt = Ee iuZt + t 0 iue -(A(t)-A(s)) M t,u s dB s .
Then, by Fubini's theorem for stochastic integrals, we easily obtain

t 0 e iuZs ds = t 0 Ee iuZs ds + t 0 dB s t s iue -(A(r)-A(s)) M r,u s dr. (3.3)
As Z t is Gaussian with variance K(0, t), it converges in distribution to a Gaussian variable of law γ = N (0, K(∞)). Because of Cesàro's result, we have

lim t→∞ 1 t t 0 Ee iuZs ds = e -u 2 2 K(∞) .
We wish to find an asymptotic equivalent to the stochastic factor of (3.3).

To this aim, consider

N u s,t (v) := t s iue A(v)-A(r) M r,u v dr. First, on the set { ∞ 0 < N u •,t (s) > s ds < ∞}, the local martingale t 0 N u s,t (
s)dB s converges a.s. to a finite variable and so, is of the order of o(t). Actually, we decompose it as

t 0 N u s,∞ (s)dB s - t 0 N u t,∞ (s)dB s . (3.4) On the set { ∞ 0 < N u •,t (
s) > s ds = ∞}, the LLN for martingales implies a.s.

t 0 dB s ∞ s iue -(A(r)-A(s)) M r,u s dr = o t 0 ∞ s iue -(A(r)-A(s)) M r,u s dr 2 ds .
Indeed, we obtain the following upper bound by using the initial definition of M r,u s :

|N u s,t (s)| ≤ |u| t s e A(s)-A(r) dr = |u|e A(s) (I t -I s )
where

I t := t 0 e -A(r) dr = I ∞ -e -A(t) a(t) + o e -A(t) a(t)
, we find by the triangle inequality that t 0 e 2A(s) (I t -I s ) 2 ds = O(t). So we have

E t 0 N u t,∞ (s)dB s 2 = t 0 E(N u t,∞ (s) 2 )ds ≤ |u| 2 t 0 e 2A(s) ds(I ∞ -I t ) 2 = O(1).
Borel-Cantelli's lemma permits to conclude that 1 t t 0 N u t,∞ (s)dB s converges a.s. to 0. Theorem 3.4. Suppose that g (t)/g 2 (t) converges to 0. Then, with probability 1, the empirical measure µ t = r r+t µ + 1 r+t t 0 δ Xs ds converges weakly to µ ∞ . Moreover, the mean µ t = 1 r+t t 0 X s ds + r r+t µ also converges almost surely.

Proof. We remind that reader that g (t)/g 2 (t) converges to 0. We start by proving that μt converges a.s. Decompose the process

µ t = µ 1 t + µ 2 t + µ 3 t
where

µ 1 t = µ + r(x -µ)H(t), µ 2 t = (H(t) -H(∞)) t 0 (r + s)e cG(s) dB s , µ 3 t = t 0 (r + s)e cG(s) (H(∞) -H(s)) dB s .
The convergence of H(t) obviously implies the convergence of µ 1 t . The deterministic factor of µ 2 t is equivalent to 1 cg(t)t 2 e -cG(t) and, due to Lemma 2.3, the quadratic variation of the stochastic factor in µ 2 t is of the order of

t 2 e 2cG(t) g(t)
. By Lemma 2.3 and the law of the iterated logarithm ( [START_REF] Lépingle | Sur le comportement aymptotique des martingales locales[END_REF] Theorem 3), we finally have µ 2 t a.s.

--→ t→∞ 0. Finally, µ 3 t is a L 2 -bounded martingale and thus converges a.s. Putting all the pieces together, we conclude that µ t a.s.

--→ t→∞ µ + H(∞)r(x -µ) + ∞ 0 (r + s)e cG(s) (H(∞) -H(s)) dB s .
To show that µ t converges, we first point out that the deterministic factor of X t converges. Decompose the process X into three parts:

X t = µ ∞ + φ(t)U t + o(1)
where

µ ∞ := x + cr(µ -x)F (∞) + ∞ 0 1 -(r + s)ce cG(s) (F (∞) -F (s)) dB s , U t := e -cG(t) r + t t 0 (r + s)e cG(s) dB s , φ(t) := c(r + t) (F (∞) -F (t)) e cG(t) .
Again, we prove the result for the Fourier transform. We have the following:

1 t t 0 e iuXs ds = e iu(µ ∞ +o(1)) t t 0 e iuφ(s)Us ds.
By Proposition 3.3, φ(t)U t is ergodic. So 1 t t 0 e iuφ(s)Us ds converges a.s. Corollary 3.5. Suppose that g (t)/g 2 (t) converges to 0 and g(∞) < ∞.

Then the limit lim µ t = µ ∞ is the Gaussian measure µ ∞ = N µ ∞ , 1 2g(∞)c . 3.3.
Asymptotic behaviour of X. We prove here that, depending on g, X exhibits three different behaviours: X converges either almost surely, or in probability (and not a.s.), or it diverges. First, we describe roughly the asymptotic behaviour of X. Proposition 3.6. Suppose that g(∞) < ∞. Then we have

P lim sup t→∞ X t = +∞ = P lim inf t→∞ X t = -∞ = 1.
Proof. Let A be a non negligible subset of R. We have the asymptotic equivalence We conclude that P(lim sup t→∞ X t = +∞) = 1. The proof is the same for lim inf t→∞ X t . Proposition 3.7. Suppose that g (t)/g 2 (t) converges to 0 and g(∞) = ∞. Then X t converges in probability to a random variable X ∞ and a.s. µ t converges weakly to δ X∞ .

Proof. As Y is Gaussian and E(Y 2 t ) = O(g(t) -1 ), it converges in L 2 and so in probability to 0. Decomposing X as X t = Y t + t 0 Y s ds r+s , it remains to show that the previous integral converges in probability. Using the explicit form of Y , Fubini's theorem for stochastic integrals ensures

t 0 Y s ds r + s = r(x-μ) t 0 e -cG(s) ds (r + s) 2 + t 0 (r+u)e cG(u) t u e -cG(s) ds (r + s) 2 dB u .
The quadratic variation of the Brownian integral converges by Lemma 2.3 and thus X converges to X ∞ in L 2 . Remark that the law of the iterated logarithm does not imply here that X converges a.s. since we do not know whether log G(t)/g(t) converges to 0 or not. We then easily have that µ t converges toward δ X∞ in probability. By Theorem 3.4, a.s. µ t converges (weakly) and so we conclude.

Proposition 3.8. Suppose that g (t)/g 2 (t) converges to 0 and g(t) -1 log G(t) is bounded for t large enough. Then there exists C > 0 such that

P(lim sup t→∞ |Y t | ≤ C) = 1.
Proof. We write Y as a Brownian (local) martingale: Y t = e -cG(t) r+t Y 0 + t 0 (r + s)e cG(s) dB s . To estimate the quadratic variation of Y , we use Lemma 2.3 and thus, by the law of the iterated logarithm, there exists C such that a.s. lim sup

t→∞ |Y t | ≤ C.
Corollary 3.9. Suppose that g(t) -1 log G(t) is lower bounded away from 0 and upper bounded for t large enough. Then X t is bounded a.s., converges in probability (but not a.s.) to X ∞ = µ ∞ and a.s. µ t converges weakly to δ X∞ .

Proof. Y is a.s. bounded and µ t converges a.s., so X t = Y t + µ t is also a.s. bounded. As Y t is Gaussian, it converges (in law) to a centered Gaussian variable. The latter being bounded, Y t converges in probability to 0. By the law of the iterated logarithm, Y t does not converge a.s. to 0 (since log G(t)/g(t) > 0 for large t). So, X t converges in probability to µ ∞ . We conclude by uniqueness of the limit that a.s. µ t converges weakly to δ µ ∞ . Proposition 3.10. Suppose that g (t)/g 2 (t) converges to 0 and lim t→∞ g(t) -1 log G(t) = 0. Then the process Y t := X t -µ t converges to 0 a.s. Moreover, both X t and µ t converge to µ ∞ a.s. and a.s. µ t converges weakly to δ µ ∞ .

Proof. We only need to prove that Y t := X t -µ t converges a.s. to 0. We have already seen that Y t = e -cG(t) r+t t 0 (r + s)e cG(s) dB s + r(x -µ) e -cG(t) r+t =: U t + v t . The deterministic term v t converges obviously to 0 and the law of the iterated logarithm implies that U t converges a.s. to 0. By uniqueness of the limit of µ t , we conclude that µ ∞ = δ µ ∞ .

Study of Y with respect to the critical points of V

From now on, we assume that g (t)/g 2 (t) converges to 0 (this hypothesis is only needed to study the behaviour of Y near a local minimum of V ). We study the process Y t = X t -μt , which is the solution to

dY t = dB t -g(t)∇V (Y t ) + Y t r + t dt; Y 0 = x -µ. (4.1)
More precisely, we study the behaviour of Y near the critical points of V . We show in particular, for each local minimum of V , that Y stays close to it with positive probability; whereas this probability is zero for any unstable critical point.

4.1.

Behaviour near the critical points of V .

Proposition 4.1. Almost surely, ∀ε > 0, ∀t > 0,

T ε t := inf{s ≥ t; d(Y s , M in ∪ M ax) < ε} < ∞. Proof. Let ε > 0. Applying Itô's formula to x → V (x), we obtain dV (Y t ) = (∇V (Y t ), dB t ) -g(t)|∇V (Y t )| 2 + 1 r + t (Y t , ∇V (Y t )) - 1 2 ∆V (Y t ) dt.
It then follows from the growth condition (2.1) that on the set {z; d(z, M in∪ M ax) > ε} and for t ≥ 0, the function y → Thus Y s∧T ε t gets close to M in ∪ M ax and there is a contradiction. Finally, P(T ε t < +∞) = 1 for all t > t 0 . For t ≤ t 0 , we conclude since t → T ε t is increasing.

, V (Y s∧T ε t ) + C(s ∧ T ε t ), s ≥ t and V (Y s∧T ε t ) + C(s ∧ T ε t ) + s∧T ε t 0 g(u) - 1 2 g(∞) |∇V (Y u∧T ε t )| 2 du

4.2.

Case of a stable critical point: local minimum. We will prove that if Y 0 is near a local minimum m, then the set {Y s ; s ≥ 0} stays near m with a positive probability. Indeed, a second-order Taylor expansion permits to compare (y -m, ∇V (y)) with |y -m| 2 and we use a comparison theorem. Let m be a local minimum of V such that ∇ 2 V (m) > 0. By Taylor's formula, there exist a > 0 and ε 0 > 0 such that for all |y -m| ≤ ε 0 we have (y -m, ∇V (y)) ≥ a|y -m| 2 . Without any loss of generality, we suppose m = 0 in the proofs.

Proposition 4.3. Suppose that g(t) -1 log G(t) is bounded on R + . Let ε 0 > ε > 0.
Then, there exists a positive stopping time T 0 such that for all T > T 0 , we have, on the event {|Y T -m| < ε}, that P (∀s ≥ 0; |Y s+T -m| < ε) > 0. Moreover, for any T > T 0 , we have on the event {∀s ≥ T ; |Y s -m| < ε}:

|Y t+T -m| = O g(t + T ) -1 log G(t + T ) a.s.
Proof. Consider the time-shifted process Y t := Y t+T . Let ε > 0. We will construct a one-dimensional process U such that for all t ≥ 0, we have a.s.

| Y t | ≤ U t .
1) Suppose that d = 1. As V (m) > 0, there exists a > 0 such that for all |y| ≤ ε, yV (y) ≥ ay 2 . Introduce the nonnegative process U , unique solution to the SDE

dU t = sign ( Y t ) dB T t -ag(t + T )U t dt + dL t , U 0 = | Y 0 |, (4.3)
where L is the local time of U in 0. Let α(t) be the function such that α(t) 0 e 2aG(s+T ) ds = t and α(0) = 0. Then, the process

A t := α(t) 0 e aG(s+T ) dL s is the local time in zero of W t = α(t) 0
e aG(s+T ) sign ( Y s ) dB T s . Denote by W + the reflected Brownian motion associated to W . Skorokhod's lemma (see [START_REF] Harrison | On the skew Brownian motion[END_REF]) then entails that e aG(α(t)+T ) U α(t) = W + t . So, the (strong) solution to (4.3) is

U t = U 0 + e -aG(t+T ) W + α -1 (t) . (4.4)
By a martingale comparison theorem, we prove that | Y t | ≤ U t a.s. (Indeed, let l be a function of class C 2 such that ∀x > 0 : l(x) > 0 and l (x) > 0, and ∀x ≤ 0 : l(x) = 0. We apply Itô's formula to l(| Y t | -U t ) to show that, on the event {| Y s | > U s }, we have l(| Y t | -U t ) ≤ 0 a.s.). Finally, as α -1 (t) = t 0 e 2aG(s+T ) ds, we conclude by the law of the iterated logarithm (LIL), that a.s.

U t = O g(t + T ) -1 log G(t + T ) . 2) Suppose that d ≥ 2. Define τ := inf{t > 0; Y t = 0}. Itô's formula implies d| Y t∧τ | = dW t∧τ -g(t ∧ τ + T ) Y t∧τ | Y t∧τ | , ∇V ( Y t∧τ ) dt - | Y t∧τ | r + t ∧ τ + T dt + d -1 2| Y t∧τ | dt where W t = t 0 Ys | Ys| , dB T s is a standard Brownian motion. The condition ∇ 2 V (0) > 0 implies
U t = O g(t + T ) -1 log G(t + T
) . Now, we prove that P (∀s ≥ T ; |Y s -m| < ε) > 0. Let τ T := inf{s > T ; |Y s -m| > ε}. For all T < t < τ T , we have a.s. |Y t -m| ≤ U t . By the LIL applied to U , we conclude that, for T large enough, P sup s≥T U s < ε > 0 and finally P(τ T = ∞) > 0.

Corollary 4.4. Suppose that g(t) -1 log G(t) converges to 0 when t tends to infinity. Then, there exists T 0 > 0 such that for all T > T 0 , the process Y t converges almost surely to m on the event {∀s ≥ T ; |Y s -m| < ε} .

4.3.

Case of an unstable critical point: local maximum or saddle point. If M is a local maximum of V, then as ∆V (M ) < 0, ε 1 := sup{ε; ∀|y| < ε, ∆V (M + y) < 0} exists and is finite. If M is a saddle point of V , then as ∇ 2 V admits a negative eigenvalue in M , there exists an unstable direction e associated to M . Let P e : R d → Re be the projection on Re. The amount ε 2 := sup{ε; ∀|y| < ε, ∂ 2 ee V (M + y) < 0 and (∂ e V (P e (y)), ∂ e V (y)) > 0} exists and is finite.

Proposition 4.5. Let M be an unstable critical point of V . If M is a local maximum, suppose that 0 < ε < ε 1 . If M is a saddle point, suppose that 0 < ε < ε 2 .
Let T be a positive stopping time such that |Y T -M | < ε. Then

P (∀s ≥ T ; |Y s -M | < ε) = 0.
Proof. Note that T < ∞ a.s. by Proposition 4.1. Suppose that M is a local maximum and M = 0, because the method of the proof is the same for M = 0 (in that case, we have an additional term M log(t + T )). Let D(t, Y t ) be the drift term of V (Y t ). On the event A := {∀s ≥ T ; |Y s | < ε}, we obtain

D(t + T, Y t+T ) = g(t + T )|∇V (Y t+T )| 2 + (Y t+T , ∇V (Y t+T )) r + t + T - 1 2 ∆V (Y t+T ) ≥ C 1 t + T + C 2
where C 1 = 1 2 inf{(y, ∇V (y)); |y| < ε} and C 2 = -1 2 sup{∆V (y); |y| < ε} > 0. We thus find for t large enough that D(t + T, Y t+T ) ≥ C > 0 and so

E (V (Y t+T )1l A ) ≤ E (V (Y T )1l A ) -CtP(A) + o(t). (4.7)
Finally, this last inequality is impossible (since V is positive) unless P(A) = 0.

Suppose now that M is a saddle point. We apply Itô's formula to x → V (P e (x)) and follow the previous computation with C 1 = 1 2 inf{(P e (y), ∂ e V (P e (y))); |y| < ε} and C 2 = -1 2 sup{∂ 2 ee V (P e (y)); |y| < ε} > 0.

Asymptotics

Through this section, we always suppose that g(∞) = +∞ and g (t)/g 2 (t) converges to 0, even if we do not remind the reader in the statements of the results. In particular, it implies that for all T > 0, G -1 (t + T ) -G -1 (t) goes to 0 when t tends to infinity. Proof. We show a stronger result: EV (Y t ) is bounded. For all n ∈ N, define the stopping time τ n = inf{t ≥ 0; |Y t | > n}. Then, there exists C > 0 such that we have by localization:

EV (Y t∧τn ) ≤ EV (Y 0 ) + e Ct < ∞.
Let n go to infinity and use Fatou's lemma to find, for all t ≥ 0, that V (Y t ) ∈ L 1 . For t large enough, we have that -g(t)V (x) + aV (x) ≤ -1 2 g(t)V (x). So, as W is c-strictly convex and by the growth hypothesis (2.1), the following holds for t large enough

d dt EV (Y t ) ≤ - 1 2 g(t)EV (Y t ). Now, solving u = -1 2 g(t)u leads to EV (Y t ) = O(1)
. In order to obtain the ergodic result for Y , we introduce a dynamical system φ, whose asymptotics are close to Y (see [START_REF] Benaïm | Dynamics of stochastic approximation algorithms[END_REF] for more details): Definition 5.2. The process Y is an asymptotic pseudotrajectory for the flow φ if

∀T > 0, ∀α > 0, lim t→+∞ P sup 0≤h≤T |Y t+h -φ h (Y t )| ≥ α = 0. ( 5
d dt φ t (x) = -∇V (φ t (x)); φ 0 (x) = x. (5.2)
Then (Y G -1 (t) , t ≥ 0) is an asymptotic pseudotrajectory for φ.

Proof.

Let Y t = Y G -1 (t) and B t = B G -1 (t)
. We will use Markov's inequality and then prove that

lim t→∞ E sup 0≤h≤T | Y t+h -φ h ( Y t )| = 0. Define κ(t) := (r + G -1 (t))g(G -1 (t)). A simple computation yields to Y t+h -φ h ( Y t ) = B t+h -B t + h 0 ∇V (φ s ( Y t )) -∇V ( Y t+s ) ds - h 0 Y t+s ds κ(t + s)
.

Applying Itô's formula to h → e -2 Ch | Y t+h -φ h ( Y t )| 2 , we have: 1 2 e 2 Ch d(e -2 Ch | Y t+h -φ h ( Y t )| 2 ) = Y t+h -φ h ( Y t ), d B t+h + Y t+h -φ h ( Y t ), Y t+h κ(t + h) dh- -C| Y t+h -φ h ( Y t )| 2 dh+ Y t+h -φ h ( Y t ), ∇V (φ h ( Y t )) -∇V ( Y t+h ) dh+ e 2 Ch 2g(G -1 (t + h))
dh.

First, we notice that G -1 (t) = 1/g(G -1 (t)). By the convexity assumption on V , we also remark that

-C| Y t+h -φ h ( Y t )| 2 + Y t+h -φ h ( Y t ), ∇V (φ h ( Y t )) -∇V ( Y t+h ) ≤ 0,
and so we deduce the following upper bound for all 0 ≤ h ≤ T :

1 2 | Y t+h -φ h ( Y t )| 2 ≤ e 2 Ch h 0 e -2 Cs ( Y t+s -φ s ( Y t ), d B t+s ) + e 2 CT 2 G -1 (t + T ) -G -1 (t) +e 2 Ch h 0 Y t+s -φ s ( Y t ), Y t+s ds κ(t + s) . (5.3)
To conclude, we will now find an upper bound for each right-hand term of (5.3). By BDG's inequality for the local martingale h 0 e -2 Cs ( Y t+sφ s ( Y t ), d B t+s ) and a rough upper bound for its quadratic variation, there exists α > 0 such that:

E sup 0≤h≤T h 0 e -2 Cs ( Y t+s -φ s ( Y t ), d B t+s ) ≤ α(G -1 (t+T )-G -1 (t))E sup 0≤h≤T | Y t+h -φ h ( Y t )| 2 1 2
.

We now estimate the remaining term of (5.3) by the triangle inequality. As κ is nondecreasing, we have the following bound by Lemma 5.1:

E T 0 | Y t+s | 2 κ(t + s) + | Y t+s -φ s ( Y t )| 2 κ(t + s) ds ≤ M T κ(t) + T κ(t) E sup 0≤h≤T | Y t+h -φ h ( Y t )| 2 .
So, we obtain for t large enough:

E sup 0≤h≤T | Y t+h -φ h ( Y t )| 2 ≤ 4e 4 CT (G -1 (t + T ) -G -1 (t)) + 4M e 2 CT T κ(t) ,
and the result follows since G -1 (t + T ) -G -1 (t) and 1/κ(t) converge to 0.

Lemma 5.4. Suppose that for all T > 0, G -1 (t+T )-G -1 (t) vanishes when t tends to infinity. Then, µ G -1

t := 1 t t 0 δ Y G -1 (s) ds, t ≥ 0 is a tight family of measures. Proof. It is enough to show that a.s. ϕ(t) := t 0 V (Y G -1 (s) )ds = O(t). Let A > 0 and K be a compact set such that ∀x ∈ K c , V (x) ≥ A. Then µ G -1 t (V ) ≥ Aµ G -1 t (1l K c ).
From the growth assumption (2.1), there exists a > 0 and for all ε > 0, there exists k ε such that ∆V ≤ aV and

V ≤ k ε + ε|∇V | 2 . It then yields ϕ(t) ≤ k ε t + ε t 0 |∇V (Y G -1 (s) )| 2 ds and t 0 ∆V (Y G -1 (s) )ds ≤ aϕ(t). (5.4) Applying Itô's formula to V (Y G -1 (t) ), we obtain V (Y G -1 (t) ) -V (Y 0 ) = G -1 (t) 0 (∇V (Y s ), dB s ) - t 0 (Y G -1 (s) , ∇V (Y G -1 (s) )) (r + G -1 (s))g(G -1 (s)) ds - t 0 |∇V (Y G -1 (s) )| 2 ds + 1 2 t 0 ∆V (Y G -1 (s) ) ds g(G -1 (s)) . (5.5)
Consider the (local-)martingale term of (5.5). On the set { ∞ 0 |∇V (Y s )| 2 ds < ∞}, it is bounded in L 2 and thus converges. Whereas on the set { ∞ 0 |∇V (Y s )| 2 ds = ∞}, the strong LLN implies that, for t large enough, a.s.

G -1 (t) 0 (∇V (Y s ), dB s ) ≤ 1 2 t 0 |∇V (Y G -1 (s) )| 2 ds.
By (5.5), we then find for t large enough

t 0 |∇V (Y G -1 (s) )| 2 ds ≤ t 0 ∆V (Y G -1 (s) )ds -2V (Y G -1 (t) ) + 2V (Y 0 ) -2 t 0 (Y G -1 (s) , ∇V (Y G -1 (s) )) (r + G -1 (s))g(G -1 (s)) ds ≤ aϕ(t) g(G -1 (t)) + 2V (Y 0 ) + O t 0 ds G -1 (s)g(G -1 (s))
.

So, we have a.s.

t 0 |∇V (Y G -1 (s) )| 2 ds = O(t) + aϕ(t).
Putting this result in (5.4) and choosing ε small enough, we conclude that ϕ(t) = O(t) a.s.

Theorem 5.5. The process Y satisfies the pointwise ergodic theorem. More precisely, there exist some (deterministic) constants a i ≥ 0, such that a i = 1 and µ Y t = 1 t t 0 δ Ys ds converges (for the weak convergence of measures) toward 1≤i≤n a i δ m i .

Proof. By Benaïm & Schreiber [START_REF] Benaïm | Weak asymptotic pseudotrajectories for semiflows: ergodic properties[END_REF], Proposition 5.3 implies that the limit points of the empirical measure of Y G -1 (t) are included in the set of all the "invariant measures" for d dt φ t (x) = -∇V (φ t (x)) with the initial condition φ 0 (x) = x. All these invariant measures are included in Span{δ

m 1 , • • • , δ mn , δ M 1 , • • • , δ Mp }. Let µ G -1 t = 1 t t 0 δ Y G -1 (
s) ds. By Lemma 5.4, the empirical measure µ G -1 t converges. One also shows that µ t is a Cauchy sequence in L 1 : there exists C > 0 such that for any s > 0,

|Eµ t+s -Eµ t | ≤ s t(t + s) t 0 E|X u |du + 1 t + s t+s t E|X u |du ≤ C s t + s . So, the limit-measure of µ G -1 t writes n i=1 a i δ m i + p i=1 b i δ M i (where a i , b i are nonnegative constants such that (a i + b i ) = 1)
. And the same result holds for µ t . Indeed, for all continuous bounded function ψ and t > s, we have (by an integration by parts)

t s ψ(Y u )du = G(t) g(t) µ G -1 G(t) ψ - G(s) g(s) µ G -1 G(s) ψ + t s g (u)G(u) g 2 (u) µ G -1 G(u) ψdu = (t -s)µ G -1 G(t) ψ + G(s) g(s) µ G -1 G(t) ψ -µ G -1 G(s) ψ + t s g (u)G(u) g 2 (u) (µ G -1 G(u) ψ -µ G -1 G(t) ψ)du.
As µ G -1 G(t) ψ converges a.s., we deduce that

µ t ψ = o(1) + µ G -1 G(t) ψ + 1 t t s g (u)G(u) g 2 (u) (µ G -1 G(u) ψ -µ G -1 G(t) ψ)du.
So, by Lemma 2.3, µ t converges. We also wish to show that b i = 0 for all i. Proposition 4.5 implies that, for an unstable critical point M , there exists a direction j such that for all ε > 0, P(∀s ≥ T, |Y

s -M (j) | ≤ ε) = 0. Consider a continuous function f , supported by a small ball (of radius α > 0) around M : f vanishes in all critical points except M and f (M ) = 1. Then, we have a.s. At this stage, we have proved that Y satisfies the pointwise ergodic theorem. The main question is whether X also satisfies the pointwise ergodic theorem or not. e i(u,Ys) e i(u,µ s ) -e i(u,µ ∞ ) ds.

The first right member converges a.s. to e i(u,µ ∞ ) 1≤p≤n e i(u,mp) . For the second right member, Cesàro asserts that it converges a.s. to 0 if and only if µ t converges a.s. 5.2. Almost sure convergence toward the local minima of V . Let 0 < ε < ε 0 and T > T 0 be as in Section 4. Let m be a local minimum of V such that |Y T -m| < ε. Lemma 5.7. If lim t→∞ g(t) -1 log G(t) = 0, then for all α > 0, we have ∞ 0 e -αg(t) dt < +∞. Proof. For all ε > 0, there exists t such that, for all s ≥ t, we have g(s) ≥ ε -1 log G(s). Moreover, there exists a > 0 such that for t large enough g(t) ≥ a and then G(t) ≥ at. So, we conclude that ∞ 1 e -ε -1 α log(at) dt < ∞ (choose for instance ε = α/2). Proposition 5.8. If g(t) -1 log G(t) converges to 0, then Y t converges a.s. and for all i, we have P lim t→∞ Y t = m i > 0 and P lim t→∞ Y t = M i = 0.

Proof. Benaïm ([1] Proposition 4.6) asserts that if -∇V (x) is a continuous globally integrable vector field, and if for all α > 0, we have ∞ 0 e -αg•G -1 (t) dt < +∞ and P(sup t |Y t | < ∞) = 1, then Y is a.s. an asymptotic pseudotrajectory for the flow induced by -∇V . Actually, the first and last conditions are fulfilled here. Moreover, as G -1 is nondecreasing, the (finite) integral ∞ 0 e -αg(t) dt is a upper bound for the preceding one. Thus, Y is a.s. an asymptotic pseudotrajectory for the flow φ defined by (5.2) and φ restricted to the limit points of Y does not admit any other attractor than the set of limit points. Finally, Y converges a.s. and its limit points are included into the set {x; ∇V (x) = 0}.

If Y converges to Y ∞ , then due to Proposition 4.5, the limit Y ∞ is not a local maximum. On the event {∀s ≥ T ; |Y s -m i | < ε}, occuring with a positive probability by Proposition 4. log G(t) = 1 a.s., we conclude that U t a.s.

--→ t→∞ 0.

Corollary 5.9. Suppose that lim g(t) -1 log G(t) = 0.

1) If V is a strictly uniformly convex function (m its unique minimum), then Y t a.s.

--→ t→∞ m.

2) A necessary condition for the almost sure convergence of Y to 0 is that the potential V admits a unique minimum at 0 (e.g. V is symmetric and strictly convex). (2) If V admits 0 as a (non-unique) local minimum, then X t converges on the event {Y t a.s.

--→ 0} and diverges elsewhere. More precisely, one has Proof. Denote m = (m (1) , • • • , m (d) ). First, suppose that m = 0. By Proposition 5.8, Y t converges toward 0 with a positive probability. On this event, Proposition 4.3 implies that the integral t 0 Ys r+s ds converges. So, µ t converges toward this (limit) integral and the result follows for X. On the other hand, if m = 0, then P(Y t → m) > 0 and so the j th -coordinate of µ t converges to sgn(m (j) )∞. So, the direction j is unstable and X t does not converge a.s. Moreover, on the set {Y ∞ = 0}, we have

µ t log t -Y ∞ ≤ 1 log t t 0 |Y s -Y ∞ | ds r + s ≤ 1 log t t 0 log G(s)
g(s) ds r + s .

The latter upper bound tends to 0 by the law of the iterated logarithm (Proposition 4.3). As Xt log t = µ t log t + Yt log t , the result follows. Remark 5.11. Any polynomial h satisfies the required condition. In particular, one can choose g(t) = t α (log(1 + t)) β with α > 0, or α = 0 and β > 2.

t 0 δ

 0 Xs (A)ds ∼ t→∞ tl where l is a positive constant depending on A. So, ∞ 0 δ Xs (A)ds = ∞ a.s. and µ ∞ is diffusive. It then implies that for all K > 0, ∞ 0 δ Xs ([K, ∞[)ds = ∞ a.s. and so P

Corollary 4 . 2 .

 42 Almost surely, the sequence of stopping times T n := inf{s > n; d(Y s , M in ∪ M ax) < ε} satisfies: T n → ∞, and ∀n ≥ 1, P(T n < +∞) = 1 and d(Y Tn , M in ∪ M ax) < ε.

5. 1 .

 1 Ergodicity. Lemma 5.1. The process Y is bounded in L 2 .

. 1 ) 5 . 3 .

 153 Proposition Let φ : R + × R d → R d be the flow generated by

t0

  1l {|Y (j) s -M (j) |≤α} ds = o(t) and 1 t t 0 f (Y s )ds converges almost surely to b. So, we conclude that b = 0.

  3, we have a.s. |Y t+T -m i | ≤ U t . As lim t→∞

5. 3 . 1 )

 31 Back to X. Theorem 5.10. Assume that g(t) -1 log G(t) = O(h(t) -2 ), with h : R + → R + such that ∞ 0 ds (1+s)h(s) < +∞. Then, on the set {Y ∞ = 0}, Xt log t converges to Y ∞ . Moreover:(If 0 is the unique local minimum of V , then

( 3 )

 3 If 0 is not a local minimum of V , thenP lim t→∞ |X t | = ∞ = P lim t→∞ |µ t | = ∞ = 1.

  1. Indeed, what we need here is not only the convergence of Y t to zero, but the convergence of the integral The process Y converges almost surely to Y ∞ , where Y ∞ belongs to the set of the local minima of V . For each local minimum m of V , one has P(Y ∞ = m) > 0. (2) On the set {Y ∞ = 0}, both X t and µ t converge almost surely to µ ∞ := µ +

	t 0 Y s	ds r+s , which depends on the speed of convergence of Y t . The main result
	of this paper is the following description of the asymptotic behaviour of X,
	shown in §5.3:
	Theorem 2.7. Suppose that g(t) -1 log G(t) = O(h(t) -1 ), where G is a
	primitive of g and	∞ 0	ds (1+s) h(s) < ∞.
	(1) ∞ 0 Y s	ds r+s , whereas on the set {Y ∞ = 0}, we have that
		lim t→∞	Xt log t = Y ∞ .

  that there exists a > 0 such that ∀|y| ≤ ε, (y, ∇V (y)) ≥ a|y| 2 .

				(4.5)
	Let us introduce the (d -1)-dimensional Bessel process R. Consider the
	time-shifted process U t := e -aG(t+T ) R t 0 e 2aG(s+T ) ds , which is the nonnegative
	(strong) solution to
	dU t = dβ T t -ag(t + T )U t dt +	d -1 2U t	dt,	(4.6)

where β is a Brownian motion. On the event {∀s ≥ T ; |Y s | < ε}, we apply the previous comparison theorem to obtain a.s. | Y t | ≤ U t . On the other hand, R is the radial part of a d-dimensional Brownian motion. So, the LIL implies a.s. R t = O( (t + T ) log log(t + T )) and

  To answer it, we remind that µ t converges a.s. if and only if Proposition 5.6. The measure µ t converges weakly if and only if µ t converges a.s.Proof. We have shown in Theorem 5.5 that Y is pointwise ergodic. Consider the Fourier transform of µ t and recall that X t = Y t + µ t . We have for all u ∈ R d :

	1 t	0	t	e i(u,Xs) ds =	e i(u,µ ∞ ) t	0	t	e i(u,Ys) ds +	1 t	0	t
	t 0 Y s a.s.	ds r+s converges, and if Y t	a.s. --→ t→∞	0 polynomially fast then µ t converges

Acknowledgements: The authors thank O. Raimond and M. Benaïm for careful reading and discussions. We are very grateful to two anonymous referees for their helpful comments. A.K. has been partially supported by the Swiss National Science Foundation grant 200020-112316/1.


Xueying Shen

Filippo Focacci
email: ffocacci@decisionbrain.com

Fabio Furini
email: fabio.furini@dauphine.fr

Virginie Gabrel
email: virginie.gabrel@dauphine.fr

Daniel Godard
email: daniel.godard@decisionbrain.com

Mathematical models for real-world production planning problems with sequence-dependent set-up costs

Introduction

Production Planning (PP) aims at making production decisions optimizing the Manufacturing and Supply Chain performances.

Given a limited amount of manufacturing resources, raw materials, and inventory space, PP minimizes the overall manufacturing and supply chain costs while meeting service level requirements and obeying operational constraints. In this research we focus on the sequencedependent set-up costs which occur when we change a product production during manufacturing. In real-world manufacturing, this cost typically depends on both the predecessor and current production. The set-up costs are also called clean-up costs, changeover costs or start-up costs (in the literature). These costs are one important part of the production costs for many real-world production systems. The set-up costs optimization may dramatically affect the quality of the production plans in terms of machine operational efficiency, cycle time, service level and inventory turns. Finally, different levels of details are typically taken into consideration which largely affect the complexity of the resulting models. In this work, we study a real-world PP problem with sequence-dependent set-up costs and we propose different formulations with different levels of details of the operational constraints.

Modelling

One of the most general problem in PP is the Capacitated Multi Product Lot-sizing Problem (CMPLP) (we refer the interested reader to [START_REF] Karimi | The capacitated lot sizing problem : a review of models and algorithms[END_REF][START_REF] Drexl | Lot sizing and scheduling : Survey and extensions[END_REF][START_REF] Pochet | Production Planning by Mixed Integer Programming[END_REF] for further details). We are given a set I = {1, 2, . . . , n} of products that must be produced during a predefined time window. For each product i and time slot t, a set-up cost sc i t , a production cost pc i t and a holding cost hc i t are prescribed. A maximum capacity p i t is given for each period and product. The problem asks for satisfying a foreseen demand d i t (per period and product) at the lowest possible cost. The formal mixed integer programming model is shown below :

min i t (pc i t x i t + sc i t y i t + hc i t s i t ) s i t-1 + x i t = d i t + s i t ∀t, ∀i x i t ≤ p i t y i t ∀t, ∀i x i t ≥ 0 ∀t, ∀i y i t ∈ {0, 1} ∀t, ∀i s i t ≥ 0 ∀t, ∀i,
where the x i t , y i t , s i t are decision variables representing the producing quantity, the decision of producing and the storage quantity of product i at time slot t respectively. In many literature works, the set-up costs depend on the current job and time slot under consideration. In this work we consider the Sequence-Dependent Set-Up Costs Production Planning Problem which requires planning of lot sizing and sequencing of production lots (see e.g. [START_REF] Haase | Capacitated lot-sizing with sequence dependent setup costs[END_REF][START_REF] Guimarãesa | Modeling lotsizing and scheduling problems with sequence dependent setups[END_REF]).

Contribution

In this paper, we present different methods to model sequence-dependent set-up costs and analyze the resulting effect. First, different ways to represent set-up costs are explored, which are shown to be memory and computation efficient in comparison to a naive matrix format. The naive matrix format considers storing the changeover costs between every pair of product productions in a matrix. Second, for production planning with time buckets, various ways of modelling and approximating sequence dependent set-up costs are provided. The models are evaluated and compared through the effect on lot sizing, production grouping and production wheel, as well as from the perspective of complexity and accuracy. More important, we propose a framework to help determining the best model to be used depending on the characteristics of the production planning problem and the time horizon. Finally, we are taking into consideration an application which comes from a real world planning project delivered by DecisionBrain. Extensive computational experiments are reported and discussed.