Mean Reflected Stochastic Differential Equations with Jumps
Philippe Briand Abir Ghannoum Céline Labart

Mean reflected SDEs with jumps]Mean Reflected Stochastic Differential
  Equations with jumps



Ph. Briand]Philippe Briand
Univ. Grenoble Alpes, Univ. Savoie Mont Blanc, CNRS, LAMA, 73000 Chambéry, France
philippe.briand@univ-smb.fr

A. Ghannoum]Abir Ghannoum
Univ. Grenoble Alpes, Univ. Savoie Mont Blanc, CNRS, LAMA, 73000 Chambéry, France
Univ. Libanaise, LaMA-Liban, P.O. Box 37, Tripoli, Liban
abir.ghannoum@univ-smb.fr

C. Labart]Céline Labart
Univ. Grenoble Alpes, Univ. Savoie Mont Blanc, CNRS, LAMA, 73000 Chambéry, France
celine.labart@univ-smb.fr





This paper is devoted to the study of reflected Stochastic Differential
Equations with jumps when the constraint is not on the paths of the solution
but acts on the law of the solution. This type of reflected equations have
been introduced recently by Briand, Elie and Hu <cit.> in
the context of BSDEs, when no jumps occur. In <cit.>, the
authors study a numerical scheme based on particle systems to approximate
these reflected SDEs. In this paper, we prove existence and uniqueness of
solutions to this kind of reflected SDEs with jumps and we generalize
the results obtained in <cit.> to this context.



[
    [
    March 23, 2018
==================






§ INTRODUCTION



Reflected stochastic differential equations have been introduced in the
pionneering work of Skorokhod (see <cit.>), and their numerical approximations by Euler
schemes have been widely studied (see <cit.>, <cit.>,
<cit.>, <cit.>, <cit.>). Reflected stochastic differential
equations driven by a Lévy process have also been studied in
the literature (see <cit.>, <cit.>). More recently, reflected
backward stochastic differential equations with jumps have been introduced and studied
(see <cit.>, <cit.>, <cit.>, <cit.>, <cit.>, <cit.>), as well as their
numerical approximation (see <cit.> and <cit.>). The main
particularity of our work comes from the fact that the constraint acts on the law of
the process X rather than on its paths. The study of such equations is
linked to the mean field games theory, which has been introduced by Lasry and
Lions (see <cit.>, <cit.>, <cit.>, <cit.>) and whose probabilistic point of view is studied in <cit.>
and <cit.>. Stochastic differential equations
with mean reflection have been
introduced by Briand, Elie and Hu in their backward forms in <cit.>. In
that work, they show that mean reflected stochastic processes exist
and are uniquely defined by the associated system of equations of the following form:

    X_t  =X_0+∫_0^t b(X_s) ds + ∫_0^t σ(X_s) dB_s  + K_t,   t≥ 0, 
       [h(X_t)] ≥ 0,   ∫_0^t [h(X_s)]   dK_s = 0,    t≥ 0.

  Due to the fact that the reflection process K depends on the law of the
  position, the authors of <cit.>, inspired by mean field games, study the convergence of a numerical scheme based
  on particle systems to compute numerically solutions to (<ref>).


  In
this paper, we extend previous results to the case of jumps, i.e. we study existence and
uniqueness of solutions to the following mean reflected stochastic
differential equation (MR-SDE in the sequel)


    X_t  =X_0+∫_0^t b(X_s^-) ds + ∫_0^t σ(X_s^-) dB_s  + ∫_0^t∫_E F(X_s^-,z) Ñ(ds,dz) + K_t,   t≥ 0, 
       [h(X_t)] ≥ 0,   ∫_0^t [h(X_s)]   dK_s = 0,    t≥ 0,

where E=ℝ^*, Ñ is a compensated Poisson measure
Ñ(ds,dz)=N(ds,dz)-λ(dz)ds, and B is a Brownian process
independent of N. We also propose a numerical scheme
based on a particle system to compute numerically solutions to (<ref>)
and study the rate of convergence of this scheme.



Our main motivation for studying (<ref>) comes from financial problems
submitted to risk measure constraints. Given any position X, its risk
measure ρ(X) can be seen as the amount of own fund needed by the investor
to hold the position. For example, we can consider the following risk measure:
ρ(X) = inf{m: 𝔼[u(m+X)]≥ p} where u is a utility
function (concave and increasing) and p is a given threshold (we refer the
reader to <cit.> and to <cit.> for more details on risk measures). Suppose that we are given a portfolio X of assets whose dynamic, when
there is no constraint, follows the jump diffusion model

    d X_t = b(X_t) d t + σ(X_t) d B_t +∫_E F(X_t-,z) Ñ(dt,dz),      t≥ 0.

 Given a risk measure ρ, one can
ask that X_t remains an acceptable position at each time t. The constraint rewrites [h(X_t)] ≥ 0 for t≥ 0 where h=u-p.

In order to satisfy this constraint, the agent has to add some cash in the portfolio through the time and the dynamic of the wealth of the portfolio becomes 

    d X_t = b(X_t) d t + σ(X_t) d B_t +∫_E F(X_t-,z) Ñ(dt,dz)+d K_t,      t≥ 0,

where K_t is the amount of cash added up to time t in the portfolio to balance the "risk" associated to X_t. Of course, the agent wants to cover the risk in a minimal way, adding cash only when needed: this leads to the Skorokhod condition [h(X_t)] d K_t = 0. Putting together all conditions, we end up with a dynamic of the form (<ref>) for the portfolio.



The paper is organized as follows. In Section <ref>,
we show that, under Lipschitz assumptions on b, σ and F and
bi-Lipchitz assumptions on h, the system admits a unique strong solution,
i.e. there exists a unique pair of process (X,K) satisfying system
(<ref>) almost surely, the process K being an increasing and
deterministic process. Then, we show that, by adding some regularity on the
function h, the Stieltjes measure dK is absolutely continuous with respect
to the Lebesgue measure and we obtain the explicit expression of its density. In Section <ref> we show that the system (<ref>) can be seen
as the limit of an interacting particles system with oblique reflection of
mean field type. This result allows to define in Section <ref> an
algorithm based on this interacting particle system together with a classical
Euler scheme which gives a strong approximation of the solution of
(<ref>). When h is bi-Lipschitz, this leads to an approximation error in L^2-sense proportional to n^-1
+ N^-1/2, where n is the number of points of
the discretization grid and N is the number of particles. When h is
smooth, we get an approximation error proportional to n^-1
+ N^-1. By the way, we improve the speed of convergence obtained in
<cit.>. Finally, we illustrate these results numerically in Section <ref>.











§ EXISTENCE, UNIQUENESS AND PROPERTIES OF THE SOLUTION.

Throughout this paper, we consider the following set of assumptions.
(A.1)
    
     0.5mm 
    (i) Lipschitz assumption: There exists a constant C_p>0, such that for all x,x'∈ℝ and p>0, we have
    
    | b(x)-b(x')|^p+|σ(x)-σ(x')|^p+∫_E| F(x,z)-F(x',z)|^pλ(dz)≤ C_p| x-x'|^p.

    (ii) The random variable X_0 is square integrable independent of B_t and N_t.


(A.2)
    
     0.5mm 
    (i) The function h:⟶ is an increasing function and there exist 0<m≤ M such that 
    ∀ x∈, ∀ y∈, m|x-y|≤|h(x)-h(y)|≤ M|x-y|.

    (ii) The initial condition X_0 satisfies: [h(X_0)]≥0.


(A.3)
    
    ∃ p>4 such that X_0 belongs to ^p: [|X_0|^p]<∞.


(A.4)
    
    The mapping h is a twice continuously differentiable function with bounded derivatives.




 §.§ Preliminary results


Define the function

    H    : ℝ×𝒫(ℝ) ∋  (x,ν) ↦∫ h(x+z) ν(dz),

and the inverse function in space of H evaluated at 0, namely:

    G̅_0    : 𝒫(ℝ) ∋ν↦inf{x  ∈ℝ : H(x,ν) ≥ 0 } ,

as well as G_0, the positive part of G̅_0:

    G_0    : 𝒫(ℝ) ∋ν↦inf{x  ≥ 0 : H(x,ν) ≥ 0 }.


We start by studying some properties of H and G_0.


Under <ref>, we have:


  (i) For all ν in 𝒫(ℝ), the mapping H(·,ν): ℝ∋ x ↦ H(x,ν) is a bi-Lipschitz function, namely: 
    
    
    ∀ x,y ∈ℝ, m|x-y| ≤ |H(x,ν)-H(y,ν)|≤ M|x-y|.

    

  (ii)  For all x in ℝ, the mapping H(x,·) : 𝒫(ℝ)∋ν↦ H(x,ν) satisfies the following Lipschitz estimate:
    
    ∀ν,ν' ∈𝒫(ℝ),  |H(x,ν)-H(x,ν')|≤|∫ h(x+·) (dν-dν')|.

    



The proof is straightforward from the definition of H (see (<ref>)).


Note that thanks to Monge-Kantorovitch Theorem, assertion (<ref>) implies that for all x in ℝ, the function H(x,·) is Lipschitz continuous w.r.t. the Wasserstein-1 distance. Indeed, for two probability measures ν and ν', the Wasserstein-1 distance between ν and ν' is defined by:

    W_1(ν,ν')=sup_φ 1-Lipschitz|∫φ(dν-dν')|=inf_X∼ν ; Y∼ν'[|X-Y|].

Therefore

    ∀ν,ν'∈𝒫(ℝ), |H(x,ν)-H(x,ν')|≤ MW_1(ν,ν').

Then, we have the following result about the regularity of G_0:


Under <ref>, the mapping G_0 :𝒫(ℝ) ∋ν↦G_0(ν) is Lipschitz-continuous in the following sense:

    |G_0(ν)-G_0( ν') | ≤1m|∫ h(G̅_0(ν)+·) (dν-dν')|,

where G̅_0(ν) is the inverse of H(·,ν) at point 0. In particular

    |G_0(ν)-G_0( ν') | ≤MmW_1(ν, ν').



The proof is given in (<cit.>, Lemma 2.5).




 §.§ Existence and uniqueness of the solution of (<ref>)


We emphasize that existence and uniqueness results hold only under <ref> which is the standard assumption for SDEs and <ref> which is the assumption used in <cit.>. The convergence of particles systems requires only an additional integrability assumption on the initial condition, namely <ref>. We sometimes add the smoothness assumption <ref> on h in order to improve some of the results.

We first recall the existence and uniqueness result of  in the case of SDEs.


A couple of processes (X,K) is said to be a flat deterministic solution to (<ref>) if (X,K) satisfy (<ref>) with K being a non-decreasing deterministic function with K_0=0.

Given this definition we have the following result.



Under Assumptions <ref> and <ref>, the mean reflected SDE (<ref>) has a unique deterministic flat solution (X, K). Moreover,

    ∀ t≥ 0, K_t=sup_s≤ tinf{x≥0:[h(x+U_s)] ≥ 0}=sup_s≤ tG_0(μ_s),
 
where (U_t)_0≤ t≤ T is the process defined by:

    U_t=X_0+∫_0^t b(X_s^-) ds + ∫_0^t σ(X_s^-) dB_s  + ∫_0^t∫_E F(X_s^-,z) Ñ(ds,dz),

and (μ_t)_0≤ t≤ T is the family of marginal laws of (U_t)_0≤ t≤ T.



The proof for the case of continuous backward SDEs is given in <cit.>. For the ease of the reader, we sketch the proof for the forward case with jumps.

Let X̂ be a given process such that, for all t>0, [sup_s≤ t|X̂_s|^2]<∞.
We set 
    Û_t=X_0+∫_0^t b(X̂_s^-) ds + ∫_0^t σ(X̂_s^-) dB_s  + ∫_0^t∫_E F(X̂_s^-,z) Ñ(ds,dz),
 and define the function K by setting 

    K_t=sup_s≤ tinf{x≥0:[h(x+Û_s)] ≥ 0}=sup_s≤ tG_0(μ̂_s).

The function K being given, let us define the process X by the formula

    X_t=X_0+∫_0^t b(X̂_s^-) ds + ∫_0^t σ(X̂_s^-) dB_s  + ∫_0^t∫_E F(X̂_s^-,z) Ñ(ds,dz)+K_t.

Let us check that (X, K) is the solution to (<ref>). By definition of K, [h(X_t)] ≥ 0 and we have dK almost everywhere,  
    K_t=sup_s≤ tinf{x≥0:[h(x+Û_s)] ≥ 0}>0,

so that [h(X_t)]=[h(Û_t+K_t)]=0   dK-a.e. since h is continuous and nondecreasing.



Next, we consider the set 𝒞^2={X , (sup_t≤ T|X_s|^2)<∞} and the map Ξ:𝒞^2⟶𝒞^2 which associates to X̂ the process X. Let us show that Ξ is a contraction. Let X̂, X̂'∈𝒞^2 be given, and define K and K' as above, using the same Brownian motion. We have from Assumption <ref>, Cauchy-Schwartz and Doob inequality


    [sup_t≤ T|X_t-X'_t|^2]    ≤ 4 [sup_t≤ T{|∫_0^t (b(X̂_s^-)-b(X̂'_s^-)) ds|^2 + |∫_0^t (σ(X̂_s^-)-σ(X̂'_s^-)) dB_s|^2 
            + |∫_0^t∫_E (F(X̂_s^-,z)-F(X̂'_s^-,z)) Ñ(ds.dz)|^2 + |K_t-K'_t|^2}]
       ≤  4 {[sup_t≤ Tt∫_0^t |b(X̂_s^-)-b(X̂'_s^-)|^2 ds] + [sup_t≤ T|∫_0^t (σ(X̂_s^-)-σ(X̂'_s^-)) dB_s|^2] 
            + [sup_t≤ T|∫_0^t∫_E (F(X̂_s^-,z)-F(X̂'_s^-,z)) Ñ(ds,dz)|^2] + sup_t≤ T|K_t-K'_t|^2}
       ≤ C {T[∫_0^T |b(X̂_s^-)-b(X̂'_s^-)|^2 ds]+ [∫_0^T |σ(X̂_s^-)-σ(X̂'_s^-)|^2ds] 
            + ∫_0^T∫_E[|F(X̂_s^-,z)-F(X̂'_s^-,z)|^2] λ(dz)ds + sup_t≤ T|K_t-K'_t|^2}
       ≤ C{T^2C_1[sup_t≤ T|X̂_t^--X̂'_t^-|^2]+ TC_1[sup_t≤ T|X̂_t^--X̂'_t^-|^2] 
            +  TC_1[sup_t≤ T|X̂_t^--X̂'_t^-|^2] + sup_t≤ T|K_t-K'_t|^2}
       ≤ C(T^2C_1+TC_2)[sup_t≤ T|X̂_t-X̂'_t|^2]+Csup_t≤ T|K_t-K'_t|^2.

   
From the representation (<ref>) of the process K and Lemma <ref>, we have that  
 

    sup_t≤ T|K_t-K'_t|^2    ≤M/m[sup_t≤ T|Û_t-Û'_t|^2]
       ≤ C(T^2C_1+TC_2)[sup_t≤ T|X̂_t-X̂'_t|^2].
  
 
Therefore,


    [sup_t≤ T|X_t-X'_t|^2]    ≤  C(1+T)T[sup_t≤ T|X̂_t-X̂'_t|^2].
   
 Hence, there exists a positive 𝒯, depending on b, σ, F and h only, such that for all T <𝒯, the map Ξ is a contraction. We first deduce the existence and uniqueness of solution on [0, 𝒯] and then on ^+ by iterating the construction.  




 §.§ Regularity results on K, X and U




Note that from this construction, we deduce that for all 0≤ s< t: 
    
    K_t-K_s
       = sup_s≤ r≤ tinf{x≥ 0 : 𝔼[ h ( x+X_s+∫_s^r b(X_u^-) du +∫_s^rσ(X_u^-) dB_u + ∫_s^r∫_E F(X_u^-,z) Ñ(du,dz) )] }.
       



From the representation (<ref>) of the process K, we have

    K_t   =sup_r≤ tinf{x≥0:[h(x+U_r)] ≥ 0}
       
    =sup_r≤ t G_0(U_r)
       
    =max{sup_r≤ sG_0(U_r),sup_s≤ r≤ t G_0(U_r)}
       
    =max{K_s,sup_s≤ r≤ t G_0(U_r)}
       
    =max{K_s,sup_s≤ r≤ t G_0(X_s-K_s+U_r-U_s)}
       
    =max{K_s,sup_s≤ r≤ t[G̅_̅0̅(X_s-K_s+U_r-U_s)^+]}.

By the definition of G̅_̅0̅, we can observe that for all y∈ℝ, G̅_̅0̅(X+y)=G̅_̅0̅(X)-y, so we get

    K_t   =max{K_s,sup_s≤ r≤ t[(K_s+G̅_̅0̅(X_s+U_r-U_s))^+]}
       
    =K_s+max{0,sup_s≤ r≤ t[(K_s+G̅_̅0̅(X_s+U_r-U_s))^+-K_s]}.

Note that sup_r (f(r)^+)=(sup_r f(r))^+=max(0,sup_r f(r)) for all function f, and obviously

    K_t   =K_s+sup_s≤ r≤ t[{(K_s+G̅_̅0̅(X_s+U_r-U_s))^+-K_s}^+]
       
    =K_s+sup_s≤ r≤ t[(G̅_̅0̅(X_s+U_r-U_s))^+]
       
    =K_s+sup_s≤ r≤ tG_0(X_s+U_r-U_s),

and so

    K_t-K_s=sup_s≤ r≤ tG_0(X_s+U_r-U_s).


In the following, we make an intensive use of this representation formula of the process K.


Let (F_s)_s≥0 be a filtration on (Ω,F,ℙ) such that (X_s)_s≥0 is an (F_s)_s≥0-adapted process.


Suppose that Assumptions <ref> and <ref> hold. Then, for all p≥ 2, there exists a positive constant K_p, depending on T, b, σ, F and h such that


  (i) [sup_t≤ T|X_t|^p]≤ K_p(1+[|X_0|^p]).

  (ii)∀ 0≤ s≤ t≤ T,    [sup_s≤ u≤ t|X_u|^p|F_s]≤ C(1+|X_s|^p).






Under the same conditions, we conclude that 
    [sup_t≤ T|U_t|^p]≤ K_p(1+[|X_0|^p]).




We have

    [sup_t≤ T|X_t|^p]   ≤ 5^p-1{|X_0|^p+sup_t≤ T(∫_0^t |b(X_s^-)| ds)^p + sup_t≤ T|∫_0^t σ(X_s^-) dB_s|^p  
         + sup_t≤ T|∫_0^t∫_E F(X_s^-,z) Ñ(ds,dz)|^p + K_T^p}.

Let us first consider the last term K_T =sup_t≤ TG_0(μ_s). From the Lipschitz property of Lemma <ref> of G_0 and the definition of the Wasserstein metric we have

    ∀ t≥ 0, |G_0(μ_t)|≤M/m[|U_t-U_0|],

since G_0(μ_0)=0 as [h(X_0)]≥ 0 and where U is defined by (<ref>). Therefore

    |K_T|^p=|sup_t≤ TG_0(μ_t)|^p   ≤ 3^p-1(M/m)^p{sup_t≤ T(∫_0^t |b(X_s^-)| ds)^p + sup_t≤ T|∫_0^t σ(X_s^-) dB_s|^p  
         + sup_t≤ T|∫_0^t∫_E F(X_s^-,z) Ñ(ds,dz)|^p},

and so

    [sup_t≤ T|X_t|^p]   ≤ C(p,M,m)[|X_0|^p+sup_t≤ T(∫_0^t |b(X_s^-)| ds)^p + sup_t≤ T|∫_0^t σ(X_s^-) dB_s|^p  
         + sup_t≤ T|∫_0^t∫_E F(X_s^-,z) Ñ(ds,dz)|^p].

Hence, using Assumption <ref>, Cauchy-Schwartz, Doob and BDG inequality we get

    [sup_t≤ T|X_t|^p]   ≤ C{[|X_0|^p] + T^p-1[∫_0^T (1+|X_s^-|)^p ds] + C_1[∫_0^T (1+|X_s^-|)^2 ds]^p/2
         + C_2[∫_0^T (1+|X_s^-|)^p ds]}
       ≤ C_1(1+|X_0|^p)+C_2∫_0^T[sup_t≤ r|X_t|^p]dr,

and from Gronwall's Lemma, we can conclude that  for all p≥ 2, there exists a positive constant K_p, depending on T, b, σ, F and h such that

    [sup_t≤ T|X_t|^p]≤ K_p(1+[|X_0|^p]).



For the first part, we have

    X_u    =U_u+K_u
       =X_s+(U_u-U_s)+(K_u-K_s)
       =X_s+∫_s^u b(X_r^-) dr + ∫_s^u σ(X_r^-) dB_r  + ∫_s^u∫_E F(X_r^-,z) Ñ(dr,dz)
            + (K_u-K_s).

Define that _s[]=[|F_s], then we get

    _s[sup_s≤ u≤ t|X_u|^p]   ≤ 5^p-1{_s[|X_s|^p]+_s[ sup_s≤ u≤ t|∫_s^u b(X_r^-) dr|^p] + _s [sup_s≤ u≤ t|∫_s^u σ(X_r^-) dB_r|^p]  
         + _s [sup_s≤ u≤ t|∫_s^u∫_E F(X_r^-,z) Ñ(dr,dz)|^p] +|K_t-K_s|^p}
       ≤  C{|X_s|^p+T^p-1∫_s^t _s[|b(X_r^-)|^p] dr+ ∫_s^t _s[|σ(X_r^-)|^p] dr  
         + ∫_s^t∫_E _s[|F(X_r^-,z)|^p]λ(dz) dr +2|K_T|^p}
       ≤ C(T){|X_s|^p+C_1∫_s^t _s[1+|X_r^-|^p]dr+2|K_T|^p}
       ≤ C_1(1+|X_s|^p)+C_2∫_s^t _s[sup_s≤ u ≤ r
    |X_u^-|^p]dr.
 
Finally, from Gronwall's Lemma, we deduce that for all 0≤ s≤ t≤ T, there exists a constant C, depending on p, T, b, σ, F and h such that
    [sup_s≤ u≤ t|X_u|^p|F_s]≤ C(1+|X_s|^p).





Let p≥2 and let Assumptions <ref>, <ref> and <ref> hold. There exists a constant C depending on p, T, b, σ, F and h such that


  (i) ∀ 0≤ s<t≤ T,    |K_t-K_s|≤ C|t-s|^(1/2).

  (ii) ∀ 0≤ s≤ t≤ T,    [|U_t-U_s|^p]≤ C|t-s|.

  (iii) ∀ 0≤ r<s<t≤ T,    [|U_s-U_r|^p|U_t-U_s|^p]≤ C|t-r|^2.





Under the same conditions, we conclude that 
    ∀ 0≤ s≤ t≤ T,    [|X_t-X_s|^p]≤ C|t-s|.




Let us recall that 
    G̅_̅0̅(X)=inf{x∈ℝ:[h(x+ X)] ≥ 0},


    G_0(X)=(G̅_̅0̅(X))^+=inf{x≥ 0:[h(x+X)] ≥ 0},

for all process X.

From Remark <ref>, we have

    K_t-K_s=sup_s≤ r≤ tG_0(X_s+U_r-U_s).

Hence, from the previous representation of K_t-K_s, we will deduce the 1/2-Hölder property of the function t⟼ K_t. Indeed, since by definition G_0(X_s)=0, if s<t, using Lemma <ref>,

    |K_t-K_s|   =sup_s≤ r≤ tG_0(X_s+U_r-U_s)
       
    =sup_s≤ r≤ t[G_0(X_s+U_r-U_s)-G_0(X_s)]
       
    =M/msup_s≤ r≤ t[|U_r-U_s|],

and so

    |K_t-K_s|   ≤ C{[sup_s≤ r≤ t|∫_s^r b(X_u^-) du|] + ([sup_s≤ r≤ t|∫_s^r σ(X_u^-) dB_u|^2])^1/2
            + ([sup_s≤ r≤ t|∫_s^r∫_E F(X_u^-,z) Ñ(du,dz)|^2])^1/2}
       ≤ C{∫_s^t[|b(X_u^-)|]du+([∫_s^t|σ(X_u^-)|^2du])^1/2
           +([∫_s^t∫_E |F(X_u^-,z) |^2λ(dz)du])^1/2}
       ≤ C{|t-s|[1+sup_u≤ T|X_u|]+|t-s|^1/2([1+sup_u≤ T|X_u|^2])^1/2}.

Therefore, if X_0∈^p for some p≥2, it follows from Proposition <ref> that 

    |K_t-K_s|≤ C|t-s|^1/2.




    
    [|U_t-U_s|^p]   ≤ 4^p-1[(∫_s^t |b(X_r^-)| dr)^p + |∫_s^t σ(X_r^-) dB_r|^p  
         +|∫_s^t∫_E F(X_r^-,z) Ñ(dr,dz)|^p ]
       ≤ Csup_0≤ r≤ t[(∫_s^r|b(X_u^-)| du)^p + |∫_s^r σ(X_u^-) dB_u|^p  
         +|∫_s^r∫_E F(X_u^-,z) Ñ(du,dz)|^p]
       ≤ C{|t-s|^p-1[∫_s^t (1+|X_u^-|)^p du] + C_1[(∫_s^t (1+|X_u^-|)^2 du)^p/2] 
         + C_2[∫_s^t (1+|X_u^-|)^p du]}
       ≤ C_1[1+sup_t≤ T|X_t|^p]|t-s|^p+C_2[(1+sup_t≤ T|X_t|^2)^p/2]|t-s|^p/2
         +C_3[1+sup_t≤ T|X_t|^p]|t-s|.

    Finally, if X_0∈^p for some p≥2, we can conclude that there exists a constant C, depending on p, T, b, σ, F and h such that
    
    ∀ 0≤ s≤ t≤ T,    [|X_t-X_s|^p]≤ C|t-s|.



Let 0≤ r<s<t≤ T, we have

    [|U_s-U_r|^p|U_t-U_s|^p]   ≤[|U_s-U_r|^p_s[|U_t-U_s|^p]]
       ≤ C[|U_s-U_r|^p{_s[|∫_s^t b(X_s^-) ds|^p]+_s[|∫_s^t σ(X_s^-) dB_s|^p]
           +_s[|∫_s^t∫_E F(X_s^-,z) dÑ(ds,dz)|^p]}].

Then, from Burkholder-Davis-Gundy inequality, we get

    [|U_s-U_r|^p|U_t-U_s|^p]   ≤ C[|U_s-U_r|^p{_s[|∫_s^t b(X_s^-) ds|^p]+(_s[∫_s^t |σ(X_s^-)|^2 ds])^p/2
           +_s[∫_s^t∫_E |F(X_s^-,z)|^p λ(dz)ds]}]
       ≤ C[|U_s-U_r|^p{|t-s|^p(1+_s[sup_s≤ u≤ t|X_u|^p])
           +|t-s|^p/2(1+_s[sup_s≤ u≤ t|X_u|^2]^p/2)+|t-s|(1+_s[sup_s≤ u≤ t|X_u|^p])}]
       ≤ C[|U_s-U_r|^p{|t-s|(1+_s[sup_s≤ u≤ t|X_u|^p])}],

thus, from (i) and Proposition <ref>, we obtain

    [|U_s-U_r|^p|U_t-U_s|^p]   ≤ C_1|t-s|[|U_s-U_r|^p]+C_2|t-s|[|U_s-U_r|^p|X_s|^p]
       ≤ C_1|t-s||s-r|+C_2|t-s|[|U_s-U_r|^p(|X_s-X_r|^p+|X_r|^p)]
       ≤ C_1|t-r|^2+C_2|t-s|[2^p-1|U_s-U_r|^p(|U_s-U_r|^p+|K_s-K_r|^p)]
           +C_3|t-s|[|U_s-U_r|^p|X_r|^p]
       ≤ C_1|t-r|^2+C_2|t-s|[|U_s-U_r|^2p]+C_3|t-s||s-r|^p/2[|U_s-U_r|^p]
           +C_4|t-s|[|U_s-U_r|^p|X_r|^p]
       ≤ C_1|t-r|^2+C_4|t-s|[|X_r|^p_r[|U_s-U_r|^p]].

Following the proof of (ii), we can also get

    _r[|U_s-U_r|^p]≤ C|s-r|(1+_r[sup_r≤ u≤ s|X_u|^p]).

Then

    [|U_s-U_r|^p|U_t-U_s|^p]   ≤ C_1|t-r|^2+C_2|t-s||s-r|[|X_r|^p(1+_r[sup_r≤ u≤ s|X_u|^p])]
       ≤ C_1|t-r|^2+C_2|t-r|^2[|X_r|^p(1+sup_r≤ u≤ s|X_u|^p)].

Under <ref>, we conclude that

    [|U_s-U_r|^p|U_t-U_s|^p]≤ C|t-r|^2,    ∀ 0≤ r<s<t≤1.





 §.§ Density of K

Let L be the linear partial operator of second order defined by

    Lf(x):=b(x)∂/∂ xf(x)+1/2σσ^*(x)∂^2/∂ x^2f(x)+∫_E(f(x+F(x,z))-f(x)-F(x,z)f'(x))λ(dz),

for any twice continuously differentiable function f.


Suppose that Assumptions <ref>, <ref> and <ref> hold and let (X,K) denote the unique deterministic flat solution to (<ref>). Then the Stieljes measure dK is absolutely continuous with respect to the Lebesgue measure (Proposition <ref>) with density

    k:ℝ^+∋ t⟼([Lh(X_t^-)])^-/[h'(X_t^-)]1_[h(X_t)]=0.



Let us admit for the moment the following result that will be useful for our proof.


The functions t⟼[h(X_t)] and t⟼[Lh(X_t)] are continuous.




If φ is a continuous function such that, for some C≥ 0 and p≥ 1,

    ∀ x∈ℝ,    |φ(x)|≤ C(1+|x|^p),
 
then the function t⟼ [φ(X_t)] is continuous.

The proof of Lemma <ref> is given in Appendix A.1. We may now proceed to the proof of Proposition <ref>.

Let t in [0,T]. For all positive r, we have

    X_t+r   =X_t+∫_t^t+r(b(X_s^-)-∫_E F(X_s^-,z)λ(dz)) ds + ∫_t^t+rσ(X_s^-) dB_s  + ∫_t^t+r∫_E F(X_s^-,z) N(ds,dz) 
           + K_t+r-K_t.

Under <ref> and thanks to Itô's formula we get

    h(X_t+r)-h(X_t)   =∫_t^t+rb(X_s^-)h'(X_s^-) ds+∫_t^t+rσ(X_s^-)h'(X_s^-) dB_s+∫_t^t+r∫_E F(X_s^-,z)h'(X_s^-) Ñ(ds,dz)
           +∫_t^t+r h'(X_s^-) dK_s+1/2∫_t^t+rσ^2(X_s^-)h”(X_s^-) ds
           +∫_t^t+r∫_E(h(X_s^-+F(X_s^-,z))-h(X_s^-)-F(X_s^-,z)h'(X_s^-))N(ds,dz)
       
    = ∫_t^t+rb(X_s^-)h'(X_s^-) ds+∫_t^t+rσ(X_s^-)h'(X_s^-) dB_s+∫_t^t+r∫_E F(X_s^-,z)h'(X_s^-) Ñ(ds.dz)
           +∫_t^t+r h'(X_s^-) dK_s+1/2∫_t^t+rσ^2(X_s^-)h”(X_s^-) ds
           +∫_t^t+r∫_E(h(X_s^-+F(X_s^-))-h(X_s^-)-F(X_s^-)h'(X_s^-))λ(dz)ds
           +∫_t^t+r∫_E(h(X_s^-+F(X_s^-,z))-h(X_s^-)-F(X_s^-,z)h'(X_s^-))Ñ(ds,dz)
       
    = ∫_t^t+rLh(X_s^-) ds+∫_t^t+r h'(X_s^-) dK_s+∫_t^t+rσ(X_s^-)h'(X_s^-) dB_s
           +∫_t^t+r∫_E(h(X_s^-+F(X_s^-,z))-h(X_s^-))Ñ(ds,dz),

where L is given by (<ref>). Thus, we obtain

    (∫_t^t+r h'(X_s^-) dK_s)= h(X_t+r)- h(X_t)-∫_t^t+rLh(X_s^-) ds.


Suppose, at the one hand, that h(X_t)>0. Then, by the continuity of t⟼ h(X_t), we get that there exists a positive R such that for all r∈[0,R], h(X_t+r)>0. This implies in particular, from the definition of K, that dK([t,t+r])=0 for all r in [0,R].

At the second hand, suppose that h(X_t)=0, then two cases arise. Let us first assume that Lh(X_t)>0. Hence, we can find a positive R' such that for all r∈[0,R'], Lh(X_t+r)>0. We thus deduce from our Assumptions and (<ref>) that h(X_t+r)>0 for all r∈(0,R']. Therefore,
dK([t,t+r])=0 for all r∈[0,R']. Suppose next that Lh(X_t)≤0. By continuity of t⟼Lh(X_t), there exists a positive R” such that for all r∈[0,R”] it holds Lh(X_t+r)≤0. Since h(X_t+r) must be positive on this set, we could have to compensate and K_t+r is then positive for all r∈[0,R”]. Moreover, the compensation must be minimal i.e. such that h(X_t+r)=0. Equation (<ref>) becomes:

    (∫_t^t+r h'(X_s^-) dK_s)=-∫_t^t+rLh(X_s^-) ds,

on [0,R”]. Dividing both sides by r and taking the limit r⟶ 0 gives (by continuity):

    dK_t=-[Lh(X_t^-)]/[h'(X_t^-)]dt.

Thus, dK is absolutely continuous w.r.t. the Lebesgue measure with density:

    k_t=([Lh(X_t^-)])^-/[h'(X_t^-)]1_[h(X_t)]=0.




This justifies, at least under the smoothness Assumption <ref> on the constraint
function h, the non-negative hypothesis imposed on h'.



Under Assumption <ref>, and by using Lemma <ref>, we obtain the continuity of the function t⟼ h(X_t).

Under the assumptions <ref>, <ref> and <ref>, we observe that x⟼Lh(X_t) is a continuous function such that, for all x∈ℝ, there exists constants C_1, C_2 and C_3 > 0,

    |b(x)h'(x)|≤ C_1(1+|x|),


    |σ^2(x)h”(x)|≤ C_2(1+|x|^2),

and

    |∫_E(h(x+F(x,z))-h(x)-F(x,z)h'(x))λ(dz)|   ≤ C_3∫_E|F(x,z)|λ(dz)
       ≤ C_3(∫_E|F(x,z)-F(0,z)|λ(dz)+∫_E|F(0,z)|λ(dz))
       ≤ C_3∫_E|x|λ(dz)+C'_3
       ≤ C_3(1+|x|).

Finally, by using Lemma <ref>, we conclude the continuity of t⟼Lh(X_t).






§ MEAN REFLECTED SDE AS THE LIMIT OF AN INTERACTING REFLECTED PARTICLES SYSTEM.


Having in mind the notations defined in the beginning of Section 2 and especially equation (<ref>), we can write the unique solution of the SDE (<ref>) as:

    X_t  =X_0+∫_0^t b(X_s^-) ds + ∫_0^t σ(X_s^-) dB_s  + ∫_0^t∫_E F(X_s^-,z) Ñ(ds,dz)+sup_s≤ tG_0(μ_s),

where μ_t stands for the law of

    U_t  =X_0+∫_0^t b(X_s^-) ds + ∫_0^t σ(X_s^-) dB_s  + ∫_0^t∫_E F(X_s^-,z) Ñ(ds,dz).

We here are interested in the particle approximation of such a system. Our candidates are the particles, for 1≤ i≤ N,

    X^i_t  = X̅^i_0+∫_0^t b(X^i_s^-) ds + ∫_0^t σ(X^i_s^-) dB^i_s  + ∫_0^t∫_E F(X^i_s^-,z) Ñ^i(ds,dz)+sup_s≤ tG_0(μ_s^N),

where B^i are independent Brownian motions, Ñ^i are independent compensated Poisson measure, X̅^i_0 are independent copies of X_0 and μ_s^N denotes the empirical distribution at time s of the particles

    U^i_t  =X̅^i_0+∫_0^t b(X^i_s^-) ds + ∫_0^t σ(X^i_s^-) dB^i_s  + ∫_0^t∫_E F(X^i_s^-,z) Ñ^i(ds,dz),   1≤ i≤ N.

namely μ^N_s=1/N∑^N_i=1δ_U^i_s.
It is worth noticing that

    G_0(μ_s^N)=inf{x≥0:1/N∑_i=1^Nh(x+U^i_s) ≥ 0}.


In order to prove that there is indeed a propagation of chaos effect, we introduce the following independent copies of X

    X̅^i_t  =X̅^i_0+∫_0^t b(X̅^i_s^-) ds + ∫_0^t σ(X̅^i_s^-) dB^i_s  + ∫_0^t∫_E F(X̅^i_s^-,z) Ñ^i(ds,dz)+sup_s≤ tG_0(μ_s),   1≤ i≤ N,

and we couple these particles with the previous ones by choosing the same Brownian motions and the same Poisson processes.

In order to do so, we introduce the decoupled particles U̅^i, 1≤ i≤ N: 

    U̅^i_t  =X̅^i_0+∫_0^t b(X̅^i_s^-) ds + ∫_0^t σ(X̅^i_s^-) dB^i_s  + ∫_0^t∫_E F(X̅^i_s^-,z) Ñ^i(ds,dz).

Note that for instance the particles U̅^i_t are i.i.d. and let us denote by μ̅^N the empirical measure associated to this system of particles.





  (i) Under our assumptions, we have [h(X̅^i_0)]=[h(X_0)]≥ 0. However, there is no reason to have 

    1/N ∑_i=1^N h(X̅^i_0) ≥ 0,

even if N is large. As a consequence,

    G_0(μ_0^N)=inf{x≥ 0:1/N∑_i=1^Nh(x+X̅^i_0) ≥ 0}

is not necessarily equal to 0. As a byproduct, we have X^i_0 = X̅^i_0 + G_0(μ_0^N) and the non decreasing process sup_s≤ tG_0(μ_s^N) is not equal to 0 at time t=0. Written in this way, the particles defined by (<ref>) can not be interpreted as the solution of a reflected SDE. To view the particles as the solution of a reflected SDE, instead of (<ref>) one has to solve

    X^i_t  = X̅^i_0 + G_0(μ_0^N) +∫_0^t b(X^i_s^-) ds + ∫_0^t σ(X^i_s^-) dB^i_s  + ∫_0^t∫_E F(X^i_s^-,z) Ñ^i(ds,dz)+ K^N_t, 
    1/N∑_i=1^Nh(X^i_t) ≥ 0,     ∫_0^t 1/N∑_i=1^Nh(X^i_t)  dK^N_s,

with K^N non decreasing and K^N_0=0. Since, we do not use this point in the sequel, we will work with the form (<ref>).



  (ii) Following the same proof of Theorem <ref>, it is easy to demonstrate the existence and the uniqueness of a solution for the particle approximated system (<ref>).



We have the following result concerning the approximation (<ref>) by interacting particles system.



Let Assumptions <ref> and <ref> hold and T>0.


  (i) Under Assumption <ref>, there exists a constant C depending on b, σ and F such that, for each j∈{1,…,N},

    [sup_s≤ T|X^j_s-X̅^j_s|^2]≤ Cexp(C(1+M^2/m^2)(1+T^2))M^2/m^2N^-1/2.



  (ii) If Assumption <ref> is in force, then there exists a constant C depending on b, σ and F such that, for each j∈{1,…,N},

    [sup_s≤ T|X^j_s-X̅^j_s|^2]≤ Cexp(C(1+M^2/m^2)(1+T^2))1+T^2/m^2(1+[sup_s≤ T|X_T|^2])N^-1.





Let t>0. We have, for r≤ t,

    |X^j_r-X̅^j_r|   ≤|∫_0^r b(X^j_s^-)- b(X̅^j_s^-)ds| + |∫_0^r (σ(X^j_s^-)-σ(X̅^j_s^-)) dB^i_s|  
           + |∫_0^r∫_E(F(X^j_s^-,z) -F(X̅^j_s^-,z)) Ñ^j(ds,dz)| +  |sup_s≤ rG_0(μ_s^N)-sup_s≤ rG_0(μ_s)|.

Taking into account the fact that

    |sup_s≤ rG_0(μ_s^N)-sup_s≤ rG_0(μ_s)|   ≤sup_s≤ r|G_0(μ_s^N)-G_0(μ_s)| 
    ≤sup_s≤ t|G_0(μ_s^N)-G_0(μ_s)|
       ≤sup_s≤ t|G_0(μ_s^N)-G_0(μ̅^N_s)|+sup_s≤ t|G_0(μ̅_s^N)-G_0(μ_s)|,

we get the inequality

    sup_r≤ t|X^j_r-X̅^j_r|   ≤ I_1+sup_s≤ t|G_0(μ_s^N)-G_0(μ̅^N_s)|+sup_s≤ t|G_0(μ̅_s^N)-G_0(μ_s)|,

where we have set

    I_1   =∫_0^t |b(X^j_s^-)- b(X̅^j_s^-)|ds + sup_r≤ t|∫_0^r (σ(X^j_s^-)-σ(X̅^j_s^-)) dB^i_s|
           + sup_r≤ t|∫_0^r∫_E(F(X^j_s^-,z) -F(X̅^j_s^-,z)) Ñ^j(ds,dz)|.


On the one hand we have, using Assumption <ref>, Doob and Cauchy-Schwartz inequalities

    [|I_1|^2]   ≤ C{[t∫_0^t|b(X^j_s^-)- b(X̅^j_s^-) |^2ds] + [∫_0^t |σ(X^j_s^-)-σ(X̅^j_s^-)|^2 ds]
           +[∫_0^t∫_E|F(X^j_s^-,z) -F(X̅^j_s^-,z)|^2 λ(dz)ds]}
       ≤ C{tC_1∫_0^t [|X^j_s-X̅^j_s|^2] ds+ C_1∫_0^t [|X^j_s-X̅^j_s|^2] ds
           +C_1∫_0^t [|X^j_s-X̅^j_s|^2] ds}
       ≤ C(1+t)∫_0^t [|X^j_s-X̅^j_s|^2] ds.

where C depends only on b, σ and F and may change from line to line.

On the other hand, by using Lemma <ref>,

    sup_s≤ t|G_0(μ_s^N)-G_0(μ̅^N_s)|
    ≤M/msup_s≤ t1/N∑_i=1^N|U^i_s-U̅^i_s|
    ≤M/m1/N∑_i=1^Nsup_s≤ t|U^i_s-U̅^i_s|,

and Cauchy-Schwartz inequality gives, since the variables are exchangeable,

    [sup_s≤ t|G_0(μ_s^N)-G_0(μ̅^N_s)|^2]
    ≤M^2/m^21/N∑_i=1^N[sup_s≤ t|U^i_s-U̅^i_s|^2]
    = M^2/m^2[sup_s≤ t|U^j_s-U̅^j_s|^2].

Since

    U^j_s-U̅^j_s=∫_0^s (b(X^j_r^-)-b(X̅^j_r^-)) dr + ∫_0^s (σ(X^j_r^-)-σ(X̅^j_r^-)) dB^j_r  + ∫_0^s∫_E (F(X^j_r^-,z)-F(X̅^j_r^-),z) Ñ^j(dr,dz),

the same computations as we did before lead to

    [sup_s≤ t|G_0(μ_s^N)-G_0(μ̅^N_s)|^2]
    ≤ C M^2/m^2(1+t)∫_0^t [|X^j_s-X̅^j_s|^2] ds.

Hence, with the previous estimates we get, coming back to (<ref>),

    [sup_r≤ t|X^j_r-X̅^j_r|^2]   ≤ K∫_0^t [|X^j_s-X̅^j_s|^2] ds+4[sup_s≤ t|G_0(μ̅_s^N)-G_0(μ_s)|^2]
       ≤ K∫_0^t [sup_r≤ s|X^j_r-X̅^j_r|^2] ds+4[sup_s≤ t|G_0(μ̅_s^N)-G_0(μ_s)|^2],

where K=C(1+t)(1+M^2/m^2). Thanks to Gronwall's Lemma

    [sup_r≤ t|X^j_r-X̅^j_r|^2]   ≤ Ce^Kt[sup_s≤ t|G_0(μ̅_s^N)-G_0(μ_s)|^2].


By Lemma <ref> we know that

    [sup_s≤ t|G_0(μ̅_s^N)-G_0(μ_s)|^2]
    ≤1m^2[sup_s≤ t|∫ h(G̅_0(μ_s)+·) (dμ̅_s^N-dμ_s)|^2],

from which we deduce that

    [sup_r≤ t|X^j_r-X̅^j_r|^2]   ≤ Ce^Kt1m^2[sup_s≤ t|∫ h(G̅_0(μ_s)+·) (dμ̅_s^N-dμ_s)|^2].


Since the function h is, at least, a Lipschitz function, we understand that the rate of convergence follows from the convergence of empirical measure of i.i.d. diffusion processes.
The crucial point here is that we consider a uniform (in time) convergence, which may possibly damage the usual rate of convergence. We will see that however here, we are able to conserve
this optimal rate. Indeed, in full generality (i.e. if we only suppose that <ref> holds) we get that:

    1m^2[sup_s≤ t|∫ h(G̅_0(μ_s)+·) (dμ̅_s^N-dμ_s)|^2]
    ≤M^2/m^2[sup_s≤ tW_1^2(μ̅_s^N,μ_s)].

Thanks to the additional Assumption <ref> and to Remark <ref> and Proposition <ref>, we will adapt and simplify the proof of Theorem 10.2.7 of <cit.> using recent results about the control Wasserstein distance of empirical measures of i.i.d. sample to the true law by <cit.>, to obtain

    [sup_s≤ 1W_1^2(μ̅_s^N,μ_s)] ≤ CN^-1/2.

The reader can refer to the work on (<cit.>, Theorem 3.2, Proof of (i)) in order to find this result.



In the case where h is a twice continuously differentiable function with bounded derivatives (i.e. under <ref>), we succeed to take benefit from the fact that μ̅^N is an empirical measure associated to i.i.d. copies of diffusion process, in particular we can get rid of the supremum in time. In view of (<ref>), we need a sharp estimate of

    [sup_s≤ t|∫ h(G̅_0(μ_s)+·) (dμ̅_s^N-dμ_s)|^2].

Let us first observe that s⟼G̅_0(μ_s) is locally Lipschitz continuous. Indeed, since by definition H(G̅_0(μ_t),μ_t)=0, if s<t, using (<ref>),


    |G̅_0(μ_s)-G̅_0(μ_t)|   ≤1/m|H(G̅_0(μ_s),μ_t)-H(G̅_0(μ_t),μ_t)|
       =1/m|H(G̅_0(μ_s),μ_t)|,
       =1/m|[h(G̅_0(μ_s)+U_t)]|
       
    =1/m|[h(G̅_0(μ_s)+U_s+∫_s^tb(X_r^-)dr+∫_s^tσ(X_r^-)dB_r+∫_s^t∫_E F(X_r^-,z)Ñ(dr,dz))]|.


We get from Itô's formula, setting

    L̅_yf(x):=b(y)∂/∂ xf(x)+1/2σσ^*(y)∂^2/∂ x^2f(x)+∫_E(f(x+F(y,z))-f(x)-F(y,z)f'(x))λ(dz),


    h(G̅_0(μ_s)+U_t)   =h(G̅_0(μ_s)+U_s)+∫_s^tb(X_r^-)h'(G̅_0(μ_s)+U_r^-)dr+∫_s^tσ(X_r^-)h'(G̅_0(μ_s)+U_r^-)dB_r
           +∫_s^t∫_E F(X_r^-,z)h'(G̅_0(μ_s)+U_r^-)Ñ(dr,dz)+1/2∫_s^tσ^2(X_r^-)h”(G̅_0(μ_s)+U_r^-)dr
           +∫_s^t∫_E m(r,z)λ(dz)dr+∫_s^t∫_E m(r,z) Ñ(dr,dz),

with

    m(r,z)=(h(G̅_0(μ_s)+U_r^-+F(X_r^-,z))-h(G̅_0(μ_s)+U_r^-)-F(X_r^-,z)h'(G̅_0(μ_s)+U_r^-)).

Thus, we obtain

    h(G̅_0(μ_s)+U_t)   
    =h(G̅_0(μ_s)+U_s)+∫_s^tL̅_X_r^-h(G̅_0(μ_s)+U_r^-)dr+∫_s^tσ(X_r^-)h'(G̅_0(μ_s)+U_r^-)dB_r
           +∫_s^t∫_E(h(G̅_0(μ_s)+U_r^-+F(X_r^-,z))-h(G̅_0(μ_s)+U_r^-))Ñ(dr,dz),

and so

    [h(G̅_0(μ_s)+U_t)]   =[h(G̅_0(μ_s)+U_s)]+∫_s^t[L̅_X_r^-h(G̅_0(μ_s)+U_r^-)]dr
       
    =H(G̅_0(μ_s),μ_s)+∫_s^t[L̅_X_r^-h(G̅_0(μ_s)+U_r^-)]dr
       
    =∫_s^t[L̅_X_r^-h(G̅_0(μ_s)+U_r^-)]dr.

Since h has bounded derivatives and sup_s≤ T|X_s| is a square integrable random variable for each T > 0 (see Proposition <ref>), the result follows easily.

Now, we deal with (<ref>).

Let us denote by ψ the Radon-Nikodym derivative of G̅_0(μ). By definition, we have, denoting by V^i the semi-martingale s⟼G̅_0(μ_s)+U̅_s^i, since U̅^i are independent copies of U,

    R_N(s):=∫ h(G̅_0(μ_s)+·) (dμ̅_s^N-dμ_s)   
    =1/N∑_i=1^N h(G̅_0(μ_s)+U̅_s^i)-[h(G̅_0(μ_s)+U_s)]
       
    =1/N∑_i=1^N {h(G̅_0(μ_s)+U̅_s^i)-[h(G̅_0(μ_s)+U̅_s^i)]}
       
    =1/N∑_i=1^N {h(V^i_s)-[h(V_s^i)]},

It follows from Itô's formula

    h(V^i_s)   =h(V^i_0)+∫_0^sh'(V^i_r^-)dG̅_0(μ_r)+∫_0^sb(X̅^i_r^-)h'(V^i_r^-)dr+∫_0^sσ(X̅^i_r^-)h'(V^i_r^-)dB^i_r
           +∫_0^s∫_EF(X̅^i_r^-,z)h'(V^i_r^-)Ñ^i(dr,dz)+1/2∫_0^sσ^2(X̅^i_r^-)h”(V^i_r^-)dr
           +∫_0^s∫_E(h(V^i_r^-+F(X̅^i_r^-,z))-h(V^i_r^-)-F(X̅^i_r^-,z)h'(V^i_r^-))λ(dz)dr
           +∫_0^s∫_E(h(V^i_r^-+F(X̅^i_r^-,z))-h(V^i_r^-)-F(X̅^i_r^-,z)h'(V^i_r^-))Ñ^i(dr,dz)
       
    =h(V^i_0)+∫_0^sh'(V^i_r^-)ψ_r dr+∫_0^sb(X̅^i_r^-)h'(V^i_r^-)dr+1/2∫_0^sσ^2(X̅^i_r^-)h”(V^i_r^-)dr
           +∫_0^s∫_E(h(V^i_r^-+F(X̅^i_r^-,z))-h(V^i_r^-)-F(X̅^i_r^-,z)h'(V^i_r^-))λ(dz)dr
           +∫_0^sσ(X̅^i_r^-)h'(V^i_r^-)dB^i_r+∫_0^s∫_E(h(V^i_r^-+F(X̅^i_r^-,z))-h(V^i_r^-))Ñ^i(dr,dz)
       
    =h(V^i_0)+∫_0^sh'(V^i_r^-)ψ_r dr+∫_0^sL̅_X̅^i_r^-h(V^i_r^-)dr+∫_0^sh'(V^i_r^-)σ(X̅^i_r^-)dB^i_r
           +∫_0^s∫_E(h(V^i_r^-+F(X̅^i_r^-,z))-h(V^i_r^-))Ñ^i(dr,dz)
       
    =h(V^i_0)+∫_0^s{h'(V^i_r^-)ψ_r +L̅_X̅^i_r^-h(V^i_r^-)}dr+∫_0^sh'(V^i_r^-)σ(X̅^i_r^-)dB^i_r
           +∫_0^s∫_E(h(V^i_r^-+F(X̅^i_r^-,z))-h(V^i_r^-))Ñ^i(dr,dz),

Taking expectation gives

    [h(V^i_s)]   
    =[h(V^i_0)]+∫_0^s[h'(V^i_r^-)ψ_r +L̅_X̅^i_r^-h(V^i_r^-)]dr
       
    =H(G̅_0(μ_0),μ_0)+∫_0^s[h'(V^i_r^-)ψ_r +L̅_X̅^i_r^-h(V^i_r^-)]dr
       
    =0+∫_0^s[h'(V^i_r^-)ψ_r +L̅_X̅^i_r^-h(V^i_r^-)]dr.

We deduce immediately that

    R_N(s)   =1/N∑_i=1^N h(V^i_0)+1/N∑_i=1^N ∫_0^sC^i(r) dr +M_N(s)+L_N(s)
       
    =1/N∑_i=1^N h(V^i_0)+ ∫_0^s(1/N∑_i=1^N C^i(r)) dr +M_N(s)+L_N(s),

where we have set 

    C^i(r)=h'(V^i_r^-)ψ_r +L̅_X̅^i_r^-h(V^i_r^-)-[h'(V^i_r^-)ψ_r +L̅_X̅^i_r^-h(V^i_r^-)],


    M_N(s)=1/N∑_i=1^N∫_0^sh'(V^i_r^-)σ(X̅^i_r^-)dB^i_r,


    L_N(s)=1/N∑_i=1^N∫_0^s∫_E(h(V^i_r^-+F(X̅^i_r^-,z))-h(V^i_r^-))Ñ^i(dr,dz).


As a byproduct,

    sup_s≤ t|R_N(s)|   ≤|1/N∑_i=1^N h(V^i_0)|+ sup_s≤ t∫_0^s|1/N∑_i=1^N C^i(r)| dr +sup_s≤ t|M_N(s)|+sup_s≤ t|L_N(s)|
       ≤|1/N∑_i=1^N h(V^i_0)|+∫_0^t|1/N∑_i=1^N C^i(r)| dr +sup_s≤ t|M_N(s)|+sup_s≤ t|L_N(s)|.


We get, using Cauchy-Schwartz inequality, since U^i and X̅^i are i.i.d,

    [sup_s≤ t|R_N(s)|^2]   ≤ 4{𝕍[1/N∑_i=1^N h(V^i_0)]+[(∫_0^t|1/N∑_i=1^N C^i(r)| dr)^2] 
           +[sup_s≤ t|M_N(s)|^2]+[sup_s≤ t|L_N(s)|^2]}
       ≤ 4{𝕍[1/N∑_i=1^N h(V^i_0)]+t[∫_0^t|1/N∑_i=1^N C^i(r)|^2 dr] 
           +[sup_s≤ t|M_N(s)|^2]+[sup_s≤ t|L_N(s)|^2]}
       
    = 4{𝕍[1/N∑_i=1^N h(V^i_0)]+t∫_0^t𝕍(1/N∑_i=1^N C^i(r)) dr 
           +[sup_s≤ t|M_N(s)|^2]+[sup_s≤ t|L_N(s)|^2]}.
 
Thus, we get

    [sup_s≤ t|R_N(s)|^2]   ≤4/N𝕍[h(V_0)]+4t/N∫_0^t𝕍(C(r)) dr +4[sup_s≤ t|M_N(s)|^2]+4[sup_s≤ t|L_N(s)|^2]
       
    =4/N𝕍[h(V_0)]+4t/N∫_0^t𝕍(h'(V_r^-)ψ_r +L̅_X_r^-h(V_r^-)) dr 
           +4[sup_s≤ t|M_N(s)|^2]+4[sup_s≤ t|L_N(s)|^2].
 
Since M_N is a martingale with

    ⟨ M_N⟩_t=1/N^2∑_i=1^N∫_0^t(h'(V^i_r^-)σ(X̅^i_r^-))^2dr,

Doob's inequality leads to

    [sup_s≤ t|M_N(s)|^2]   ≤ 4[|M_N(t)|^2]
       =4/N^2∑_i=1^N∫_0^t[(h'(V^i_r^-)σ(X̅^i_r^-))^2]dr
       =4/N∫_0^t[(h'(V_r^-)σ(X_r^-))^2]dr.

Then, using Doob inequality for the last martingale L_N, we obtain

    [sup_s≤ t|L_N(s)|^2]   ≤ 4[|L_N(t)|^2]
       =4/N^2∑_i=1^N[(∫_0^t∫_E(h(V^i_r^-+F(X̅^i_r^-,z))-h(V^i_r^-))Ñ^i(dr,dz))^2]
           +8/N^2∑_1≤ i<j≤ N[∫_0^t∫_E(h(V^i_r^-+F(X̅^i_r^-,z))-h(V^i_r^-))Ñ^i(dr,dz)
           ∫_0^t∫_E(h(V^j_r^-+F(X̅^j_r^-,z))-h(V^j_r^-))Ñ^j(dr,dz)]
       =4/N^2∑_i=1^N∫_0^t∫_E[(h(V^i_r^-+F(X̅^i_r^-,z))-h(V^i_r^-))^2]λ(dz)dr+0
       =4/N∫_0^t∫_E[(h(V^i_r^-+F(X̅^i_r^-,z))-h(V^i_r^-))^2]λ(dz)dr.

Finally, using the fact that h has bounded derivatives, b, σ and F are Lipschitz, we get

    [sup_s≤ t|R_N(s)|^2]≤ C(1+t^2)(1+[sup_s≤ t|X_s|^2])N^-1.

This gives the result coming back to (<ref>).





§ A NUMERICAL SCHEME FOR MRSDE.


We are interested in the numerical approximation of the SDE (<ref>) on [0,T]. Here are the main steps of the scheme. Let 0=T_0<T_1<⋯<T_n=T be a subdivision of [0,T]. Given this subdivision, we denote by "_" the mapping s↦s=T_k if s∈[T_k,T_k+1), k∈{0,⋯ ,n-1}. For simplicity, we consider only the case of regular subdivisions: for a given integer n, T_k=kT/n, k=0,…,n.

Let us recall that we proved in the previous section that the particles system, for 1≤ i≤ N,

    X^i_t=X̅^i_0+∫_0^t b(X^i_s^-) ds + ∫_0^t σ(X^i_s^-) dB^i_s  + ∫_0^t∫_E F(X^i_s^-,z) Ñ^i(ds,dz) + sup_s≤ tG_0(μ^N_s),

where we have set 

    μ^N_t     =1/N∑^N_i=1δ_U^i_t, 
    
    	U^i_t     =X̅^i_0+∫_0^t b(X^i_s^-) ds + ∫_0^t σ(X^i_s^-) dB^i_s  + ∫_0^t∫_E F(X^i_s^-,z) Ñ^i(ds,dz),   1 ≤ i≤ N,

B^i being independent Brownian motions, N^i being independent Poisson processes and X̅^i_0 being independent copies of X_0, converges toward the solution to (<ref>). Thus, the numerical approximation is obtained by an Euler scheme applied to this particles system. We introduce the following discrete version of the particles system: for 1≤ i≤ N,

    X̃^i_t=X̅^i_0+∫_0^t b(X̃^i_s^-) ds + ∫_0^t σ(X̃^i_s^-) dB^i_s  + ∫_0^t∫_E F(X̃^i_s^-,z) Ñ^i(ds,dz) + sup_s≤ tG_0(μ̃^N_s),

with the notation 

    μ̃^N_t    =1/N∑^N_i=1δ_Ũ^i_t, 
    Ũ^i_t     =X̅^i_0+∫_0^t b(X̃^i_s^-) ds + ∫_0^t σ(X̃^i_s^-) dB^i_s  + ∫_0^t∫_E F(X̃^i_s^-,z) Ñ^i(ds,dz),   1 ≤ i≤ N.




 §.§ Scheme.

Using the notations given above, the result on the interacting system of mean
reflected particles of the MR-SDE of Section 3 and Remark <ref>, we deduce the following algorithm for the numerical approximation of the MR-SDE:


We emphasize that, at each step k of the algorithm, we approximate the increment of the reflection process K by the increment of the approximation:

    Δ_k K̂^N:=sup_l≤ kG_0(μ̃^N_T_l)-sup_l≤ k-1G_0(μ̃^N_T_l).



First, we consider the special case when the SDE is defined by

    X_t  =X_0+∫_0^t b(X_s^-) ds + ∫_0^t σ(X_s^-) dB_s  + ∫_0^t F(X_s^-) dÑ_s + K_t,   t≥ 0, 
       [h(X_t)] ≥ 0,   ∫_0^t [h(X_s)]   dK_s = 0,    t≥ 0.

where N is a Poisson process with intensity λ, and Ñ_t=N_t-λ t.


As suggested in Remark <ref>, the increment (<ref>) can be approached by:



Δ_k K^N:=

    inf{x≥ 0:1/N∑_i=1^N h(   x+(X̃_T_k-1^μ̃^N)^i+T/n(b((X̃_T_k-1^μ̃^N)^i)-λ F((X̃_T_k-1^μ̃^N)^j)+√(T)/√(n)σ((X̃_T_k-1^μ̃^N)^i)G^i
       +F((X̃_T_k-1^μ̃^N)^i)H^i)≥0},

where G^j∼𝒩(0,1) and H^j∼𝒫(λ(T/n)) and are i.i.d.


Indeed, using the same kind of arguments as in the sketch of the proof of Theorem <ref>, one can show that the increments of the approximated reflection process are equals to the approximation of the increments:

    ∀ k∈{1,⋯ n}:  Δ_k K^N=Δ_kK̂^N.




Returning to the general case (<ref>), we can see in <cit.>, N={N(t):=N(E×[0,t])} is a stochastic process with intensity λ that counts the number of jumps until some given time. The Poisson random measure N(dz,dt) generates a sequence of pairs {(ι_i,ξ_i),i∈{1,2,⋯,N(T)}} for a given finite positive constant T if λ<∞. Here {ι_i,i∈{1,2,⋯,N(T)}} is a sequence of increasing nonnegative random variables representing the jump times of a standard Poisson process with intensity λ, and {ξ_i,i∈{1,2,⋯,N(T)}} is a sequence of independent identically distributed random variables, where ξ_i is distributed according to f(z), where λ(dz)dt=λ f(z)dzdt.
Then, the numerical approximation can equivalently be the following form

    X̅^j_T_k   =X̅^j_T_k-1+T/n(b(X̅^j_T_k-1)-∫_Eλ F(X̅^j_T_k-1,z)f(z)dz)+√(T/n)σ(X̅^j_T_k-1)G^j
           +∑_i=H^j_T_k-1+1^H^j_T_kF(X̅^j_T_k-1,ξ_i)+Δ_k K̂^N,

Δ_kK̂^N=Δ_k K^N=

    inf{x≥ 0:1/N∑_i=1^N h(   x+X̅^j_T_k-1+T/n(b(X̅^j_T_k-1)-∫_Eλ F(X̅^j_T_k-1,z)f(z)dz)+√(T/n)σ(X̅^j_T_k-1)G^j
           +∑_i=H^j_T_k-1+1^H^j_T_kF(X̅^j_T_k-1,ξ_i))≥0},

where G^j∼𝒩(0,1) and H^j∼𝒫(λ(T/n)) and are i.i.d.



 §.§ Scheme error.






  (i) Let T>0, N and n be two non-negative integers and let Assumptions <ref>, <ref> and <ref> hold. There exists a constant C, depending on T, b, σ, F, h and X_0 but independent of N, such that: for all i=1,…,N

    [sup_s≤ t|X^i_s-X̃^i_s|^2]≤ C(n^-1+N^-1/2).


  (ii) Moreover, if Assumption <ref> is in force, there exists a constant C, depending on T, b, σ, F, h and X_0 but independent of N, such that: for all i=1,…,N

    [sup_s≤ t|X^i_s-X̃^i_s|^2]≤ C(n^-1+N^-1).





Let us fix i∈1,…,N and T>0. We have, for t≤ T,

    |X^i_s-X̃^i_s|   ≤|∫_0^s b(X^i_r^-)- b(X̃^i_r^-)dr | + |∫_0^s (σ(X^i_r^-)-σ(X̃^i_r^-)) dB^i_r|  
           + |∫_0^s∫_E(F(X^i_r^-,z) -F(X̃^i_r^-,z)) Ñ^i(dr,dz)| + sup_r≤ s|G_0(μ_r^N)-G_0(μ̃^N_r)|.

Hence, using Assumption <ref>, Cauchy-Schwartz, Doob and BDG inequalities  we get:

    [sup_s≤ t|X^i_s-X̃^i_s|^2]    ≤ 4[sup_s≤ t{|∫_0^s(b(X^i_r^-)- b(X̃^i_r^-)) dr|^2 + |∫_0^s (σ(X^i_r^-)-σ(X̃^i_r^-)) dB^i_r|^2  
           + |∫_0^s∫_E(F(X^i_r^-,z) -F(X̃^i_r^-,z)) Ñ^i(dr,dz)|^2 + sup_r≤ s|G_0(μ_r^N)-G_0(μ̃^N_r)|^2}] 
       ≤ C{[t∫_0^t|b(X^i_s^-)- b(X̃^i_s^-) |^2ds] + [∫_0^t |σ(X^i_s^-)-σ(X̃^i_s^-)|^2 ds]
           +[∫_0^t∫_E|F(X^i_s^-,z) -F(X̃^i_s^-,z)|^2 λ(dz)ds] + [sup_s≤ t|G_0(μ_s^N)-G_0(μ̃^N_s)|^2]}
       ≤ C{TC_1∫_0^t [|X^i_s-X̃^i_s|^2] ds+ C_1∫_0^t [|X^i_s-X̃^i_s|^2] ds
           +C_1∫_0^t [|X^i_s-X̃^i_s|^2] ds+ [sup_s≤ t|G_0(μ_s^N)-G_0(μ̃^N_s)|^2]}
       ≤ C∫_0^t [|X^i_s-X̃^i_s|^2] ds +4[sup_s≤ t|G_0(μ_s^N)-G_0(μ̃^N_s)|^2].

Denoting by  (μ^i_t)_0≤ t≤ T the family of marginal laws of (U^i_t)_0≤ t≤ T and (μ̃^i_t)_0≤ t≤ T the family of marginal laws of (Ũ^i_t)_0≤ t≤ T, we have

    [sup_s≤ t|G_0(μ_s^N)-G_0(μ̃^N_s)|^2]   ≤3{[sup_s≤ t|G_0(μ_s^N)-G_0(μ_s^i)|^2]
    +sup_s≤ t|G_0(μ_s^i)-G_0(μ̃^i_s)|^2
           +[sup_s≤ t|G_0(μ̃^i_s)-G_0(μ̃^N_s)|^2]},

and from Lemma <ref>,

    ≤3{1m^2[sup_s≤ t|∫ h(G̅_0(μ^i_s)+·) (dμ_s^N-dμ_s^i)|^2]
    +(M/m)^2sup_s≤ t W_1^2(μ_s^i,μ̃^i_s)
           +1m^2[sup_s≤ t|∫ h(G̅_0(μ̃^i_s)+·) (dμ̃^N_s-dμ̃^i_s)|^2]}
       ≤ C{[sup_s≤ t|∫ h(G̅_0(μ^i_s)+·) (dμ_s^N-dμ_s^i)|^2]
    +sup_s≤ t[|U_s^i-Ũ^i_s|^2]
           +[sup_s≤ t|∫ h(G̅_0(μ̃^i_s)+·) (dμ̃^N_s-dμ̃^i_s)|^2]}.


Following the Proof of (i) in Theorem <ref>, we obtain

    [sup_s≤ t|∫ h(G̅_0(μ^i_s)+·) (dμ_s^N-dμ_s^i)|^2]≤ C[sup_s≤ t W_1^2(μ_s^N,μ_s^i)]≤ CN^-1/2,


    [sup_s≤ t|∫ h(G̅_0(μ̃^i_s)+·) (dμ̃^N_s-dμ̃^i_s)|^2]≤ C[sup_s≤ t W_1^2(μ̃^i_s,μ̃^N_s)]≤ CN^-1/2.

From which, we can derive the inequality

    [sup_s≤ t|G_0(μ_s^N)-G_0(μ̃^N_s)|^2]   ≤ C_1sup_s≤ t[|U_s^i-Ũ^i_s|^2]+C_2N^-1/2
       ≤ C_1{sup_s≤ t[|U_s^i-Ũ^i_s|^2]+sup_s≤ t[|Ũ^i_s-Ũ^i_s|^2]}
    +C_2N^-1/2.

For the first term of the right hand side, we can observe that,

    sup_s≤ t[|U_s^i-Ũ^i_s|^2]]   ≤[sup_s≤ t|U_s^i-Ũ^i_s|^2]]
       ≤ 3[sup_s≤ t{|∫_0^s(b(X^i_r^-)- b(X̃^i_r^-)) dr|^2 + |∫_0^s (σ(X^i_r^-)-σ(X̃^i_r^-)) dB^i_r|^2  
         + |∫_0^s∫_E(F(X^i_r^-,z) -F(X̃^i_r^-,z)) Ñ^i(dr,dz)|^2}] 
       ≤ C{[t∫_0^t|b(X^i_s^-)- b(X̃^i_s^-) |^2ds] + [∫_0^t |σ(X^i_s^-)-σ(X̃^i_s^-)|^2 ds]
           +[∫_0^t∫_E|F(X^i_s^-,z) -F(X̃^i_s^-,z)|^2 λ(dz)ds]}
       ≤ C{TC_1∫_0^t [|X^i_s-X̃^i_s|^2] dr+ 2C_1∫_0^t [|X^i_s-X̃^i_s|^2] dr}
       ≤ C∫_0^t [|X^i_s-X̃^i_s|^2] ds.

Using Assumption <ref>, the second term sup_s≤ t[|Ũ^i_s-Ũ^i_s|^2] becomes

    sup_s≤ t[|Ũ^i_s-Ũ^i_s|^2]   ≤ 3sup_s≤ t{[|∫_s^s b(X̃^i_r^-) dr|^2 
    + |∫_s^s σ(X̃^i_r^-) dB^i_r|^2 
    + |∫_s^s ∫_E F(X̃^i_r^-,z) Ñ^i(dr,dz)|^2]}
       ≤ 3sup_s≤ t{[| b(X̃^i_s) |^2 |s-s|^2
    + |σ(X̃^i_s) |^2 |B^i_s-B^i_s|^2 
    + ∫_s^s ∫_E |F(X̃^i_r^-,z)|^2 λ(dz)dr]}
       ≤ 3sup_s≤ t{[| b(X̃^i_s) |^2 |s-s|^2
    + |σ(X̃^i_s) |^2 |B^i_s-B^i_s|^2 
    +  C∫_s^s(1+|X̃^i_r^-|^2) dr]}
       ≤ 3sup_s≤ t{(T/n)^2[|sup_s≤ r≤ s b(X̃^i_r^-) |^2]
    +  [|B^i_s-B^i_s|^2] [|σ(X̃^i_s) |^2]
    
           +C(T/n)[sup_s≤ r≤ s (1+|X̃^i_r|^2)]}
       ≤ C_1(T/n)^2[sup_s≤ T|b(X̃^i_s) |^2]
    + C_2 sup_s≤ t[|B^i_s-B^i_s|^2] [sup_s≤ T|σ(X̃^i_s) |^2]
           
    + C_3 (T/n)[sup_s≤ T (1+|X̃^i_s|^2)]
       ≤ C_1(T/n)^2(1+[sup_s≤ T|X̃^i_s |^2])
    + C_2 sup_s≤ t[|B^i_s-B^i_s|^2] (1+[sup_s≤ T|X̃^i_s |^2])
           
    + C_3  (T/n)(1+[sup_s≤ T|X̃^i_s |^2]),

and from Proposition <ref>, we get

    sup_s≤ t[|Ũ^i_s-Ũ^i_s|^2]   ≤ C_1(T/n)+ C_2 sup_s≤ t[|B^i_s-B^i_s|^2].

Then, by using BDG inequality, we obtain

    sup_s≤ t[|B^i_s-B^i_s|^2]
    =sup_s≤ t[(∫_s^s dB^i_u)^2]
    ≤sup_s≤ t|s-s|
    ≤T/n.

Therefore, we conclude 

    sup_s≤ t[|Ũ^i_s-Ũ^i_s|^2]   ≤ C_1n^-1+ C_2 n^-1
       ≤ Cn^-1,

from which we derive the inequality

    [sup_s≤ t|G_0(μ_s^N)-G_0(μ̃^N_s)|^2]   ≤ C{n^-1+N^-1/2+∫_0^t [|X^i_s-X̃^i_s|^2] ds},

and taking into account (<ref>) we get

    [sup_s≤ t|X^i_s-X̃^i_s|^2]    ≤ C{n^-1+N^-1/2+∫_0^t [|X^i_s-X̃^i_s|^2] ds}.

Since 

    [|X^i_s-X̃^i_s|^2]   ≤ 2[|X^i_s-X̃^i_s|^2]+2[|X̃^i_s-X̃^i_s|^2]
       
    =2[|X^i_s-X̃^i_s|^2]+2[|Ũ^i_s-Ũ^i_s|^2],

it follows from (<ref>) and (<ref>) that

    [sup_s≤ t|X^i_s-X̃^i_s|^2]    ≤ C{n^-1+N^-1/2+∫_0^t [|X^i_s-X̃^i_s|^2] ds}.

and finally, we conclude the proof of (i) with Gronwall's Lemma.


Following the proof of (ii) in Theorem <ref>, we obtain

    [sup_s≤ t|∫ h(G̅_0(μ^i_s)+·) (dμ_s^N-dμ_s^i)|^2]≤ CN^-1,


    [sup_s≤ t|∫ h(G̅_0(μ̃^i_s)+·) (dμ̃^N_s-dμ̃^i_s)|^2]≤ CN^-1.
 
According to the same strategy applied to proof of (i) in Theorem <ref>, the result follows easily:

    [sup_s≤ t|X^i_s-X̃^i_s|^2]≤ C(n^-1+N^-1).






Let T > 0, N and n be two non-negative integers. Let assumptions <ref>, <ref> and <ref> hold.


  (i) There exists a constant C, depending on T, b, σ, F, h and X_0 but independent of N, such that: for all i=1,…,N,

    [sup_t≤ T|X̅^i_t-X̃^i_t|^2]≤ C(n^-1+N^-1/2).


  (ii) If in addition <ref> holds, there exists a positive constant C, depending on T, b, σ, F, h and X_0 but independent of N, such that: for all i=1,…,N,

    [sup_t≤ T|X̅^i_t-X̃^i_t|^2]≤ C(n^-1+N^-1).




The proof is straightforward writing

    |X̅^i_t-X̃^i_t|≤|X̅^i_t-X^i_t|+|X^i_t-X̃^i_t|,

and using Theorem <ref> and Proposition <ref>.





§ NUMERICAL ILLUSTRATIONS.

Throughout this section, we consider, on [0,T] the following sort of processes:

    X_t  =X_0-∫_0^t (β_s+a_s X_s^-) ds + ∫_0^t (σ_s+γ_s X_s^-) dB_s  + ∫_0^t ∫_E c(z)(η_s+θ_s X_s^-) Ñ(ds,dz) + K_t,   t≥ 0, 
       [h(X_t)] ≥ 0,   ∫_0^t [h(X_s)]   dK_s = 0,    t≥ 0.
 

where (β_t)_t≥0, (a_t)_t≥0, (σ_t)_t≥0, (γ_t)_t≥0, (η_t)_t≥0 and (θ_t)_t≥0 are bounded adapted processes. This sort of processes allow us to make some explicit computations leading us to illustrate the algorithm. Our results are presented for different diffusions and functions h that are summarized below.

Linear constraint.
We first consider cases where h:ℝ∋ x ⟼ x-p∈ℝ.



  Case (i) Drifted Brownian motion and compensated Poisson process: β_t=β>0, a_t=γ_t=θ_t=0, σ_t=σ>0, η_t=η>0, X_0=x_0≥ p, c(z)=z and 
    f(z)=1/√(2π z)exp(-(ln z)^2/2)1_{0<z}.

We have 
    K_t=(p+β t-x_0)^+,

and 
    X_t=X_0-(β+λ√(e))t+σ B_t+∑_i=0^N_tηξ_i+K_t,
 
where N_t∼𝒫(λ t) and ξ_i∼ lognormal(0,1).


  Case (ii) Black and Scholes process: β_t=σ_t=η_t=0, a_t=a>0, γ_t=γ>0, θ_t=θ>0, c(z)=δ_1(z). Then


 K_t=ap(t-t^*)1_t≥ t^*, where t^*=1/a(ln(x_0)-ln(p)),

and

    X_t=Y_t+Y_t∫_0^tY_s^-1dK_s,

where Y is the process defined by:

    Y_t=X_0exp(-(a+γ^2/2+λθ)t+γ B_t)(1+θ)^N_t.



Nonlinear constraint.
Secondly, we illustrate the case of non-linear function h:

    h:ℝ∋ x⟼ x+αsin(x)-p∈ℝ,  -1<α<1,

and we illustrate this case with 


  Case (iii) Ornstein Uhlenbeck process: β_t=β>0, a_t=a>0, γ_t=θ_t=0, σ_t=σ>0, η_t=η>0, X_0=x_0 with x_0>|α|+p, c(z)=δ_1(z). We obtain

    dK_t=e^-atdsup_s≤ t(F^-1_s(0))^+,

where for all t in [0,T],

    F_t:ℝ∋ x⟼   {e^-at(x_0-β(e^at-1/a)+x)
    +αexp(-e^-atσ^2/2asinh(at))
    
       ×[1/2(exp(λ t(e^iη-1))+exp(λ t(e^-iη-1)))sin(e^-at(x_0-(β+λη)(e^-at-1/a)+x))
    
       +1/2i(exp(λ t(e^iη-1))-exp(λ t(e^-iη-1)))cos(e^-at(x_0-(β+λη)(e^-at-1/a)+x))]
    
       -p}




These examples have been chosen in such a way that we are able to give an analytic form of the reflecting process K. This enables us to compare numerically the “true” process K and its empirical approximation K̂ . When an exact simulation of the underlying process is available, we
compute the approximation rate of our algorithm.





 §.§  Proofs of the numerical illustrations.

In order to have closed, or almost closed, expression for the compensator K we introduce the process Y solution to the non-reflected SDE

    Y_t  =X_0-∫_0^t (β_s+a_s Y_s^-) ds + ∫_0^t (σ_s+γ_s Y_s^-) dB_s  + ∫_0^t∫_E c(z)(η_s+θ_s Y_s^-) Ñ(ds,dz).

By letting A_s=∫_0^ta_sds and applying Itô’s formula on e^A_tX_t and e^A_tY_t, we get

    e^A_tX_t   =X_0+∫_0^t e^A_sX_s a_s ds+∫_0^t e^A_s(-β_s-a_s X_s^-) ds + ∫_0^t e^A_s(σ_s+γ_s X_s^-) dB_s  
           + ∫_0^t∫_E e^A_s c(z)(η_s+θ_s X_s^-) Ñ(ds,dz)+∫_0^te^A_sdK_s
       
    =X_0-∫_0^t e^A_sβ_s ds + ∫_0^t e^A_s(σ_s+γ_s X_s^-) dB_s  + ∫_0^t∫_E e^A_s c(z)(η_s+θ_s X_s^-) Ñ(ds,dz)+∫_0^te^A_sdK_s,

in the same way,

    e^A_tY_t=X_0-∫_0^t e^A_sβ_s ds + ∫_0^t e^A_s(σ_s+γ_s Y_s^-) dB_s  + ∫_0^t∫_E e^A_s c(z)(η_s+θ_s Y_s^-) Ñ(ds,dz),

and so

    X_t=Y_t+e^-A_t∫_0^te^A_sdK_s+ e^-A_t∫_0^t e^A_sγ_s (X_s^-+Y_s^-) dB_s  +e^-A_t∫_0^t∫_E e^A_s c(z)θ_s(X_s^-+Y_s^-) Ñ(ds,dz).




In all of cases, we have a_t=a i.e. A_t=at, so we get

    [Y_t]   =[e^-at(x_0-∫_0^t e^asβ ds+ ∫_0^t e^as(σ_s+γ_s Y_s^-) dB_s  + ∫_0^t∫_E e^as c(z)(η_s+θ_s Y_s^-) Ñ(ds,dz))]
       
    =e^-at(x_0-∫_0^t e^asβ ds)
       
    =e^-at(x_0-β(e^at-1/a)).






From Proposition <ref> and Remark <ref>, we have

    k_t   =β1_(X_t)=p
       
    =β1_(Y_t)+K_t-p=0
       
    =β1_x_0-β t+K_t-p=0,

so, we obtain that

    K_t   =∫_0^t k_t dt
       
    =∫_0^t β1_K_t=p+β t-x_0 dt,

and as K_t≥ 0, we conclude that

    K_t=(p+β t-x_0)^+.

Next, we have 

    f(z)=1/√(2π z)exp(-(ln z)^2/2),

the density function of a lognormal random variable,
so we can obtain 

    ∫_E η z λ(dz)=λη∫_E z f(z) dz=λη(ξ)

where ξ∼ lognormal(0,1),
and we conclude that 

    ∫_E η z λ(dz)=λη√(e).

Finally, we deduce the exact solution

    X_t=X_0-(β+λ√(e)) t+σ B_t+∑_i=0^N_tηξ_i+K_t,

where N_t∼𝒫(λ t) and ξ_i∼ lognormal(0,1).



In this case, and using the same Proposition and Remark, we have

    k_t   =((-aX_t))^- 1_(X_t)=p,

which

    (X_t)=p   ⟺(Y_t)-p+e^-at∫_0^t e^asdK_s=0
       ⟺ -x_0e^-at+p=e^-at∫_0^t e^asdK_s
       ⟺ K_s=ap,

and 

    K_t≥0   ⟺ -x_0e^-at+p≥0
       ⟺ e^-at≤p/x_0
       ⟺ t≥1/a(ln(x_0)-ln(p)):=t^*.

So, we conclude that
K_t=ap(t-t^*)1_t≥ t^*, where t^*=1/a(ln(x_0)-ln(p)).

Next, by the definition of the process Y_t in this case,

    dY_t  =-a Y_t^- dt + γ Y_t^- dB_s  + θ Y_t^- dÑ_t,

we have 

    Y_t=X_0exp(-(a+γ^2/2+λθ)t+γ B_t)(1+θ)^N_t.

Thanks to Itô formula we get

    d(1/Y_t)   =-1/Y_t^2dY_t+1/2(2/Y_t^3)γ^2Y_t^2dt+d∑_s≤ t(1/Y_s^-+Δ Y_s-1/Y_s^-+1/Y_s^-^2Δ Y_s)
       
    =a/Y_tdt-γ/Y_tdB_t-θ/Y_t^-dÑ_t+γ^2/Y_tdt+d∑_s≤ t(1/(1+θ)Y_s^--1/Y_s^-+θ/Y_s^-),

and so

    dY_t^-1   
    =(a+γ^2)Y_t^-1dt-γ Y_t^-1dB_t-θ Y_t^-^-1dÑ_t+(θ^2/1+θ)d∑_s≤ tY_s^-^-1
       
    =(a+γ^2+λθ^2/1+θ)Y_t^-1dt-γ Y_t^-1dB_t-(θ/1+θ) Y_t^-^-1dÑ_t.

Then, using integration by parts formula, we obtain

    d(X_t Y_t^-1)   
    =X_t^- dY_t^-1+Y_t^-^-1dX_t+d[X,Y^-1]_t 
       
    =(a+γ^2)X_tY_t^-1dt-γ X_tY_t^-1dB_t-θ X_t^-Y_t^-^-1dÑ_t+(θ^2/1+θ)d∑_s≤ tX_s^-Y_s^-^-1
           -aX_t Y_t^-1dt+γ X_t Y_t^-1dB_t+θ X_t^- Y_t^-^-1dÑ_t+Y_t^-1dK_t
           -γ^2X_tY_t^-1dt-(θ^2/1+θ)d∑_s≤ tX_s^-Y_s^-^-1
       
    =Y_t^-1dK_t.

Finally, we deduce that 

    X_t=Y_t+Y_t∫_0^tY_s^-1dK_s.





In that case, we have

    Y_t   =e^-at(x_0-β(e^at-1/a))+ σ_s e^-at∫_0^t e^asdB_s  + e^-at∫_0^t η_s e^as  dÑ_s
    
       =e^-at(x_0-(β+λη) (e^at-1/a))+ σ_s e^-at∫_0^t e^asdB_s  + e^-at∫_0^t η_s e^as  dN_s
    
       :=f_t+G_t+F_t,

and

    X_t=Y_t+e^-atK̅_t,    K̅_t=∫_0^te^asdK_s.

Hence

    h(X_t)   =Y_t+e^-atK̅_t+αsin(Y_t+e^-atK̅_t)-p
       
    =Y_t+e^-atK̅_t+α(sin(Y_t)cos(e^-atK̅_t)+cos(Y_t)sin(e^-atK̅_t))-p
       
    =Y_t+e^-atK̅_t+α[cos(e^-atK̅_t){sin(f_t)cos(G_t)cos(F_t)+cos(f_t)sin(G_t)cos(F_t)
           +cos(f_t)cos(G_t)sin(F_t)-sin(f_t)sin(G_t)sin(F_t)}+sin(e^-atK̅_t){cos(f_t)cos(G_t)cos(F_t)
           -sin(f_t)sin(G_t)sin(F_t)-sin(f_t)cos(G_t)sin(F_t)-cos(f_t)sin(G_t)sin(F_t)}]-p.

On one side, since G_t is a centered gaussian random variable with variance V given by

    V=σ^21-e^-2at/2a=σ^2e^-atsinh(at)/a,
 
we obtain that 
    [e^iG_t]=e^-V/2,


    [sin(G_t)]=[e^iG_t-e^-iG_t/2i]=0,

and

    [cos(G_t)]=[e^iG_t+e^-iG_t/2]=(e^iG_t)=exp(-e^-atσ^2/2asinh(at))=:g(t).


On the other side, 

    [e^iF_t]=[exp(iη e^-at∫_0^t e^asdN_s)],

by taking ‘a’small, we get

    [e^iF_t]   ≈[exp(iη∫_0^t dN_s)]
       ≈[exp(iη N_t)]
       ≈exp(λ  t(e^iη-1)),

and so

    [sin(F_t)]≈exp(λ  t(e^iη-1))-exp(λ  t(e^-iη-1))/2i=:m(t),



    [cos(F_t)]≈exp(λ  t(e^iη-1))+exp(λ  t(e^-iη-1))/2=:n(t).


Using Remark <ref>, we conclude that, for small ‘a’,

    [h(X_t)]   ≈[Y_t]+e^-atK̅_t+α(g(t)m(t)cos(f_t+e^-atK̅_t)+g(t)n(t)sin(f_t+e^-atK̅_t))-p
       :=F_t(K̅_t).

 
Therefore,

    K̅_t=sup_s≤ t(F_s^-1(0))^+    and     dK_t=e^-atdsup_s≤ t(F_s^-1(0))^+.






 §.§ Illustrations.


This computation works as follows. Let 0 = T_0 < T_1 <⋯< T_n = T be
a subdivision of [0, T ] of step size T/n, n being a positive integer, let X be the unique solution of the MRSDE (<ref>) and let, for a given i, (X̃^i_T_k)_0≤ k≤ n be its numerical approximation given by Algorithm 1. For a given integer L, we draw (X̅^l)_0≤ l≤ L and (X̃^i,l)_0≤ l≤ L, L independent copies of X and X̃^i. We then approximate the 𝕃^2-error of Theorem <ref> by:

    Ê=1/L∑_l=1^Lmax_0≤ k≤ n|X̅^l_T_k-X̃^i,l_T_k|^2.





























§ APPENDICES



 §.§ Proof of Lemma <ref>

Let s and t in [0,T] such that s≤ t.

Firstly, we suppose that φ is a continuous function with compact support. In this case, there exists a sequence of Lipschitz continuous functions φ_n with compact support which converges uniformly to φ.
Therefore, by using Proposition <ref>, we get

    |[φ(X_t)]-[φ(X_s)]|   ≤ |[φ(X_t)]-[φ_n(X_t)]|+|[φ_n(X_t)]-[φ_n(X_s)]|+|[φ_n(X_s)]-[φ(X_s)]|
       ≤[|(φ-φ_n)(X_t)|]+C_n[|X_t-X_s|]+|[(φ_n-φ)(X_s)]|
       ≤ 2[∥φ_n-φ∥_∞]+C_n([|X_t-X_s|^2])^1/2
       ≤ 2[∥φ_n-φ∥_∞]+C_n|t-s|^1/2.

Thus, we obtain that

    lim sup_t→ s|[φ(X_t)]-[φ(X_s)]|≤ 2[∥φ_n-φ∥_∞].

This result is true for all n≥1, so we deduce that

    lim sup_t→ s|[φ(X_t)]-[φ(X_s)]|=0,

then we conclude the continuity of the function t⟼ [φ(X_t)].



Secondly, we consider the case that φ is a continuous function such that 
    ∀ x∈ℝ, ∃ C∈ℝ,  φ(x)≤ C(1+|x|^p).

We define a sequence of functions φ_n, such that for all n≥1 and x∈ℝ, 

    φ_n(x)=φ(x)θ_n(x)
 
with θ_n smooth such that 

    θ_n(x)=
        
          1    if |x|≤ n 
    
          0    if |x|> n+1

 Based on this definition, φ_n is a continuous function with compact support. Then we get 
 
 
    |[φ(X_t)]-[φ(X_s)]|    ≤[|φ-φ_n|(X_t)]  + |[φ_n(X_t)]-[φ_n(X_s)]|  + [|φ-φ_n|(X_s)]
    
        ≤[2|φ(X_t)|_|X_t|> n]+|[φ_n(X_t)]-[φ_n(X_s)]|+ [ 2|φ(X_s)|_|X_s|> n] 
    
        ≤ C[ (1+sup_t≤ T|X_t|^p)_sup_t≤ T|X_t|> n]+|[φ_n(X_t)]-[φ_n(X_s)]|.

 
  Thus, by using the first part of this Lemma, we obtain that
  
    lim sup_t→ s|[φ(X_t)]-[φ(X_s)]|≤ C[(1+sup_t≤ T|X_t|^p)_sup_t≤ T|X_t|> n].

  This result is true for all n≥1, then by using the dominated convergence theorem, we deduce that
  
    lim sup_t→ s|[φ(X_t)]-[φ(X_s)]|=0,

  and we conclude the continuity of the function t⟼ [φ(X_t)].






BCdRGL16

[ADEH99]ADEH99
Philippe Artzner, Freddy Delbaen, Jean-Marc Eber, and David Heath.
Coherent measures of risk.
Math. Finance, 9(3):203–228, 1999.

[BCdRGL16]BCGL17
Philippe Briand, Paul-Eric Chaudru de Raynal, Arnaud Guillin, and Céline
  Labart.
Particles systems and numerical schemes for mean reflected stochastic
  differential equations.
arXiv:1612.06886, 2016.
Submitted.

[BEH18]BEH16
Philippe Briand, Romuald Elie, and Ying Hu.
BSDEs with mean reflexion.
Ann. Appl. Probab., 28(1):482–510, 2018.

[CD18a]CD18a
R. Carmona and F. Delarue.
Probabilistic Theory of Mean Field Games with Applications I.
Springer, 2018.

[CD18b]CD18b
R. Carmona and F. Delarue.
Probabilistic Theory of Mean Field Games with Applications II.
Springer, 2018.

[CM08]CM08
Stéphane Crépey and Anis Matoussi.
Reflected and doubly reflected BSDEs with jumps: a priori estimates
  and comparison.
Ann. Appl. Probab., 18(5):2041–2069, 2008.

[DL16a]DL16a
R. Dumitrescu and C. Labart.
Numerical approximation of doubly reflected bsdes with jumps and rcll
  obstacles.
Journal of Math. Anal. and Appl., 442(1):206–243, 2016.

[DL16b]DL16b
R. Dumitrescu and C. Labart.
Reflected scheme for doubly reflected bsdes with jumps and rcll
  obstacles.
Journal of Comput. and Appl. Math., 296:827–839, 2016.

[EHO05]EHO05
E.H. Essaky, N. Harraj, and Y. Ouknine.
Backward stochastic differential equation with two reflecting
  barriers and jumps.
Stochastic Analysis and Applications, 23:921–938, 2005.

[Ess08]Ess08
E.H. Essaky.
Reflected backward stochastic differential equation with jumps and
  RCLL obstacle.
Bulletin des Sciences Mathématiques, (132):690–710, 2008.

[FG15]FG15
Nicolas Fournier and Arnaud Guillin.
On the rate of convergence in wasserstein distance of the empirical
  measure.
Probab. Theory Related Fields, (162(3-4)):707––738, 2015.

[FS02]FS02
Hans Föllmer and Alexander Schied.
Convex measures of risk and trading constraints.
Finance Stoch., 6(4):429–447, 2002.

[HH06]HH06
S. Hamadène and M. Hassani.
BSDEs with two reacting barriers driven by a Brownian motion and
  an independent Poisson noise and related Dynkin game.
Electronic Journal of Probability, 11:121–145, 2006.

[HO03]HO03
S. Hamadène and Y. Ouknine.
Reflected backward stochastic differential equation with jumps and
  random obstable.
Elec. Journ. of Prob., 8:1–20, 2003.

[KH92]KH92
Arturo Kohatsu-Higa.
Reflecting stochastic differential equations with jumps.
Technical report, 1992.

[Lep95]Lep95
D. Lepingle.
Euler scheme for reflected stochastic differential equations.
Math. Comput. Simulation, 38:119–126, 1995.

[LL06a]LL06b
J.-M. Lasry and P.-L. Lions.
Jeux à champ moyen. i. le cas stationnaire.
C.R. Math. Acad. Sci. Paris, 343(9):619–625, 2006.

[LL06b]LL06a
J.-M. Lasry and P.-L. Lions.
Jeux à champ moyen. ii. horizon fini et contrôle optimal.
C.R. Math. Acad. Sci. Paris, 343(10):679–684, 2006.

[LL07a]LL07a
J.-M. Lasry and P.-L. Lions.
Large investor trading impacts on volatility.
Ann. Inst. H. Poincaré, Anal. Non Linéaire, 24(2):311–323,
  2007.

[LL07b]LL07b
J.-M. Lasry and P.-L. Lions.
Mean field games.
Japanese Journal Math., 2:229–260, 2007.

[MR85]MR85
José-Luis Menaldi and Maurice Robin.
Reflected diffusion processes with jumps.
The Annals of Probability, 13(2):319–341, 1985.

[Pet95]Pet95
R. Petterson.
Approximations for stochastic differential equations with reflecting
  convex boundaries.
Stoc. Proc. Appl, 59(2):295–308, 1995.

[Pet97]Pet97
R. Petterson.
Penalization schemes for reflecting stochastic differential
  equations.
Bernoulli, 3(4):403–414, 1997.

[QS14]QS14
M.C. Quenez and A. Sulem.
Reflected BSDEs and robust optimal stopping for dynamic risk
  measures with jumps.
Stoc. Proc. Appl., (124(9)):3031–3054, 2014.

[Sko61]Sko61
A.V. Skorokhod.
Stochastic equations for diffusion processes in a bounded region.
Theory Probab. Appl, 6:264–274, 1961.

[Slo94]Slo94
L. Slominski.
On approximation of solutions of multidimensional sdes with
  reflecting boundary conditions.
Stoc. Proc. Appl, 50(2):197–219, 1994.

[Slo01]Slo01
L. Slominski.
Euler's approximations of solutions of sdes with reflecting boundary.
Stoc. Proc. Appl, 92(2):317–337, 2001.

[STR98]RR98
Ludger Rüschendorf S. T. Rachev.
Mass transportation problems. Vol. 1: Theory. Vol. 2:
  Applications.
1998.

[YS12]YS12
Hui Yu and Minghui Song.
Numerical solutions of stochastic differential equations driven by
  poisson random measure with non-lipschitz coefficients.
Journal of Applied Mathematics, 2012, 2012.



