Frédérique Charles

Bruno Després

Michel Mehrenberger

Enhanced convergence estimates for semi-lagrangian schemes Application to the Vlasov-Poisson equation

We prove enhanced error estimates for high order semi-lagrangian discretizations of the Vlasov-Poisson equation. It provides new insights into optimal numerical strategies for the numerical solution of this problem. The new error estimate

is based on advanced error estimates for semi-lagrangian schemes, also equal to shifted Strang's schemes, for the discretization of the advection equation. 1 p 1 ≤ p < ∞, and z L ∞ (Ω) sup x∈Ω |z(x)|.

Introduction

The aim in this work is to prove enhanced error estimates for high order semilagrangian discretizations of the Vlasov-Poisson equation [START_REF] Schaeffer | Global existence of smooth solutions to the Vlasov poisson system in three dimensions Communications in Partial Differential Equations[END_REF]. In this work enhanced means two things. First it means enhanced with respect to the time step ∆t, that is we do not want the error estimate to be spoiled by a 1 ∆t which is very often encountered in the numerical analysis of the Vlasov-Poisson equation as in [START_REF] Falcone | Optimal discretization steps in semi-Lagrangian approximation of first order PDEs[END_REF][START_REF] Pinto | Convergence of an adaptive semi-Lagrangian scheme for the Vlasov-Poisson system[END_REF][START_REF] Besse | Convergence of classes of high-order semilagrangian schemes for the Vlasov-Poisson system[END_REF]. Second, in the case of advection equation in dimension one, it means enhanced with respect to the norm, that is the error and the regularity of the solution are evaluated in the same norm.

Our main result is stated in Theorem 1. Estimate [START_REF] Falcone | Optimal discretization steps in semi-Lagrangian approximation of first order PDEs[END_REF] has the advantage that it is non singular for small ∆t. It is a strong improvement with respect to the literature [START_REF] Falcone | Optimal discretization steps in semi-Lagrangian approximation of first order PDEs[END_REF][START_REF] Pinto | Convergence of an adaptive semi-Lagrangian scheme for the Vlasov-Poisson system[END_REF][START_REF] Besse | Convergence of classes of high-order semilagrangian schemes for the Vlasov-Poisson system[END_REF]. This new error estimate is based on new estimates for the advection equation on a regular periodic grid for which we use sharp properties of semi-lagrangian schemes also equal to shifted Strang's stencil [START_REF] Strang | Trigonometric polynomials and difference methods of maximum accuracy[END_REF][START_REF] Iserles | The optimal accuracy of difference schemes[END_REF][START_REF] Iserles | Order stars[END_REF], see also [START_REF] Delpino | Arbitrary high-order schemes for linear advection and waves equations: application to hydrodynamics and aeroacoustic[END_REF][START_REF] Del Pino | 3D finite volume simulation of acoustic waves in the earth atmosphere[END_REF][START_REF] Després | Finite Volume Transport Schemes[END_REF][START_REF] Després | Uniform asymptotic stability of Strang's explicit compact schemes for linear advection[END_REF]. We will make use of the connection of these numerical schemes with B-Splines techniques [START_REF] Meinardus | On the Chebyshev Norm of Polynomial B-Splines[END_REF][START_REF] Schumaker | Spline functions: basic theory[END_REF]

with ∂E ∂x (t, x) = ˆ+∞ -∞ f (t, x, v)dv -1, (2)
where f (t, x, v) is the distribution function of charged particles (ions or electrons) and E(t, x) the self-consistent electric field. We consider periodic boundary conditions on the variable x ∈ [0, 1], that is

f (t, 0, v) = f (t, L, v), v ∈ R, t ≥ 0, (3)
and

E(t, 0) = E(t, L). ( 4
)
Assume the exact solution is sufficiently smooth and has compact support. Consider a semi-lagrangian scheme of order p direction by direction for the discretization of (63-66) with Strang's splitting in time. Assume the grid is such that ∆v = α∆x is proportional to ∆x, α kept constant. Assume n∆t ≤ T . There exists a constant C > 0 which depends on T , on the regularity of the solution and of the parameters of the problem such that the numerical error

e n = f n -f (t n ) is bounded in the L 2 norm 1 by ||e n || L 2 ≤ C min ∆x ∆t , 1 ∆x p + ∆t 2 . ( 6
)
The organization of this work is as follows. First in Section 2 we detail optimal numerical strategies which are immediate consequence of the error estimate [START_REF] Falcone | Optimal discretization steps in semi-Lagrangian approximation of first order PDEs[END_REF]. The next section is devoted to the design of enhanced error estimates for 1 We will use natural notations for the norm of discrete functions over the real line R discretized with a mesh length ∆x > 0. For example the L p norm of a discrete function w = (w i ) i∈Z defined over the entire real line is

w p = ∆x i |w i | p 1 p 1 ≤ p < ∞,
and

w ∞ = sup i |w i |.
If the domain is finite, for example Ω =]0, 1[ 2 , then N ∆x = 1 is required for some N ∈ N: in this case the discrete function is w = (w ij ) 1≤i,j≤N with norms defined by

w p = ∆x 1≤i,j≤N |w ij | p 1 p 1 ≤ p < ∞, and
w ∞ = sup 1≤i,j≤N |w ij |. (5)
These notations are compatible with the standard definition of the L p norm of a function

z L p (Ω) = ˆΩ |z(x)| p dx
the numerical discretization of the advection equation in dimension one with semi-lagrangian schemes of arbitrary orders (also equal to shifted Strang's stencils). After that we prove the main theorem. In the appendix we provide the reader with more advanced formulas for the numerical approximation of the advection equation.

Application

The computation of a numerical solution to the Vlasov-Poisson equation is very heavy task from the CPU point of view. It is therefore of major interest to study optimal numerical strategies such the upper bound ( 6) of the error is as small as possible as a function of ∆x and ∆t. We consider three different strategies with p ≥ 3 since it corresponds to the situation we are interested in. It will appear that the fact that singularity for very small ∆t is removed in the right hand side of ( 6) due to the term min ∆x ∆t , 1 has an immediate consequence on the optimal scaling laws. Essentially we will obtain that ∆t ≈ ∆x p 2 << ∆x which yields that min ∆x ∆t , 1 = 1. Note that the results would not be the same, if we would consider higher order discretizations in time (see [START_REF] Falcone | Optimal discretization steps in semi-Lagrangian approximation of first order PDEs[END_REF], for a discussion).

Minimizing the error at given storage

The storage is the memory requirement needed to run a given computation. Since the storage is proportional to 1/∆x∆v ≈ 1/∆x 2 , it means that we look for an optimal time step ∆t which minimizes the error at a given (but small) ∆x. So we are looking for a scaling law ∆t = ∆x β , β > 0 such that the right hand side (6) is as small as possible when ∆x goes to zero. Since min ∆x ∆t ,

1 ∆x p + ∆t 2 ≤ ∆x p + ∆t 2 it is immediate that β ≥ p 2 ≥ 3 2 > 1.
In this case min ∆x ∆t , 1 = min ∆x ∆x β , 1 = 1 for small ∆x. We obtain asymptotically ||e n || L 2 ≤ C∆x p with a larger constant. Our interest being to have nevertheless the biggest time step to minimize the overall cost of the computation we obtain the scaling law

β = p 2 .

Minimizing the error at given CPU cost

The CPU cost is proportional to the total number of cells and to the number of time steps, that is after normalization

CPU ≈ 1 ∆t∆x 2 .
Therefore ∆t = 1 α∆x 2 where α >> 1 is the normalized numerical value of the CPU cost. Plugging in the error formula [START_REF] Falcone | Optimal discretization steps in semi-Lagrangian approximation of first order PDEs[END_REF] we find

Error ≈ C min(α∆x 3 , 1)∆x p + 1 α 2 ∆x 4 .
The minimum is approximatively obtained by equating the two contributions

1 = min(α∆x 3 , 1) α∆x 3 2 ∆x p-2 . Therefore 1 ≤ α 2 ∆x p+4 and ∆x ≥ α -2 p+4 . So α∆x 3 ≥ α 1-6 p+4 ≥ α 1-6
3+4 >> 1 since p ≥ 3 by hypothesis. So min(α∆x 3 , 1) = 1. In conclusion the optimal strategy that minimizes the error at given CPU cost proportional to α is

∆x ≈ α -2 p+4 and ∆t ≈ α -1+ 4 p+4 = α -p p+4 .
We observe ∆t ≈ ∆x p 2 that for large p: this is the scaling law of the previous strategy.

Minimizing the CPU cost a given error

The CPU cost is a strictly convex function with respect to (∆t, ∆x). Assume the error ε is very small, that is

min ∆x ∆t , 1 ∆x p + ∆t 2 = ε, with ε << 1.
We formulate the problem as a minimum problem with a constraint. We distinguish two cases.

First case: ∆t > ∆x: Minimum solutions, if they exist, are the critical point of the Lagrangian

L = 1 ∆t∆x 2 -λ ∆x p+1 ∆t + ∆t 2 -ε
where λ is the Lagrange multiplier. The optimality conditions are

   ∂ ∆t L = -1 ∆t 2 ∆x 2 -λ -∆x p+1 ∆t 2 + 2∆t = 0, ∂ ∆x L = -2 ∆t∆x 3 -λ(p + 1) ∆x p ∆t = 0.
Since p ≥ 3 one may approximate

- ∆x p+1 ∆t 2 + 2∆t ≈ 2∆t
for small ∆x. Therefore the optimality conditions imply

-1 ∆t 2 ∆x 2 ≈ λ2∆t, -2 ∆t∆x 3 = λ(p + 1) ∆x p ∆t .
It yields

2∆t ∆x ≈ p + 1 2 ∆x p ∆t 2 that is ∆t ≈ p + 1 4 1 3 ∆x p+1 3 << ∆x since p ≥ 3.
It is in contradiction with the hypothesis ∆t < ∆x. It means that there is no solution in this case.

Second case: ∆t ≤ ∆x: The Lagrangian is

L = 1 ∆t∆x 2 -λ ∆x p + ∆t 2 -ε .
The optimality conditions are

∂ ∆t L = -1 ∆t 2 ∆x 2 -λ2∆t = 0, ∂ ∆x L = -2 ∆t∆x 3 -λp∆x p-1 = 0
whose solution is given by ∆t = p + 1 4

1 2 ∆x p 2 << ∆x.
Once again this scaling law is very close to the first one.

Exponential integrators in time

Suppose that for each splitting step, instead of solving over ∆t the 1D advection equation, one solves N times the same equation with time step ∆t/N , and that we consider the limit as N goes to infinity of such scheme. In the previous analysis [START_REF] Besse | Convergence of classes of high-order semilagrangian schemes for the Vlasov-Poisson system[END_REF], we cannot conclude that the scheme converges. On the contrary, in our new framework, one obtains the rate O(∆x p + ∆t 2 ).

A pedagogical example for the advection equation

We consider in this section the numerical discretization of the advection equation with initial condition

∂ t u + v∂ x u = 0, v > 0, u(0, x) = u 0 (x) (7)
in one dimension by means of a semi-lagrangian scheme. We recall that the solution of problem [START_REF] Hundsdorfer | Numerical solution of time-dependent advection-diffusion-reaction equations[END_REF] verifies

u(t + s, x) = u(t, x -vs) ∀t ≥ 0, s ≥ 0. ( 8
)
We introduce a Cartesian regular discretization {x i } i∈Z = {i∆x} i∈Z of R and a time step ∆t. For the sake of simplicity the numerical initial condition is always exact, that is

u 0 j = u 0 (x j ) ∀j ∈ Z. (9)
Assuming that we know the approximation of the solution of (7) at the time t n = n∆t on the mesh, we set for all i

u n+1 i = ũn i (x i -v∆t) ( 10
)
where ũn is a continuous function obtain by means of an interpolation of values {u n i } i∈{k1,...,k2} . Such a procedure is called semi-lagrangian because the discrete solution is given on a eulerian fixed mesh, but on the other hand the scheme relies on the construction of ũn (x i -v∆t) which is used in the lagrangian approximation [START_REF] Besse | Convergence of classes of high-order semilagrangian schemes for the Vlasov-Poisson system[END_REF]. This is described in Figure 1.

In order to fully understand the main difficulty about the design of optimal estimates, we consider, in a first part for pedagogical purposes, a second order lagrangian interpolation (that is k 2 = k 1 + 1). It will enlighten the interest of having an interpretation of such scheme both in terms of a semi-lagrangian interpolation and in terms of a standard numerical approximation of a partial differential equation. Then we will develop both approaches for general high order semi-lagrangian schemes.

Numerical scheme

The function ũn is obtained by a lagrangian interpolation from the two neighboring points. Let x i -v∆t be the foot of the characteristics and r ∈ Z such that

x i+r ≤ x i -v∆t < x i+r+1 , that is r ≤ -v ∆t ∆x < r + 1.
In other terms we set

r = E -v ∆t ∆x ≤ -1. (11)
Since v > 0 then r + 1 ≤ 0 for any ∆t > 0 and ∆x > 0. This is illustrated in Figure 1.

The second order reconstruction of ũn on [x i+r , x i+r+1 [ is

ũn i (x) = x i+r+1 -(x i -v∆t) ∆x u n i+r + (x i -v∆t) -x i+r ∆x u n i+r+1 . ( 12
)
Introducing this expression in [START_REF] Besse | Convergence of classes of high-order semilagrangian schemes for the Vlasov-Poisson system[END_REF] one obtains where we set

u n+1 i = νu n i+r + (1 -ν)u n i+r+1 , ( 13
) r=-2 r+1=-1 t n (n+1) ∆ x ∆ t t
ν := r + 1 + v∆t ∆x . ( 14
)
We deduce from the definitions of r (11) that

0 < ν ≤ 1. (15)
The standard way to evaluate the error in space attached to this procedure is the following. First the interpolation error of ( 12) is O(∆x 2 ) for a smooth function. Second it is without to say that ( 13) is stable in the maximum norm. So one just sums the errors from one time step to the following and obtain

max i |u(n∆t, x i ) -u n i | ≤ n O(∆x 2 ) ≤ CT ∆x 2 ∆t , n∆t ≤ T. ( 16
)
The constant C is independent of ∆t, ∆x and v. The estimate ( 16) is second order in space but is clearly non optimal for small ∆t. Indeed ( 16) corresponds to the upwind scheme for r -1: in this case it is well known that the error can be bounded by O(∆x). This analysis enlightens the fact that a naive interpolation error estimate is not optimal in our context.

The truncation error

We now desire to explain how to recast the previous analysis of the interpolation error so as to obtain optimal estimate for small ∆t. An idea is to rewrite [START_REF] Thomée | Besov spaces and applications to difference methods for initial value problems[END_REF] under the form of the Finite Difference scheme but with skewed discrete derivatives

u n+1 i -u n i+r+1 ∆t + v + (r + 1)∆x ∆t u n i+r+1 -u n i+r ∆x = 0. ( 17
)
It is usual to introduce the truncation error

ε n = (ε n i ) i∈Z ε n i = u ((n + 1)∆t, x i ) -u (n∆t, x i+r+1 ) ∆t + v + (r + 1)∆x ∆t u (n∆t, x i+r+1 ) -u (n∆t, x i+r ) ∆x , ( 18
)
where u is the solution of equation ( 7). The next task consists is showing that it can be related to a more standard truncation error.

The unit time ∆x v is the one needed for a characteristics to travel in a cell of length ∆x, from one side to the other. Let ∆s ≥ 0 be defined as ∆t minus a d ∆x v where d + 1 is the largest as possible integer. A convenient definition of ∆s is ∆s :

= (r + 1) ∆x v + ∆t. ( 19
)
We deduce from the definition of r (11) and ν (14) :

0 < ∆s ≤ ∆t, ν = v∆s ∆x . ( 20
)
Let ε be the truncation error of the same scheme but with the time step ∆s

ε n i+r+1 = u (n∆t + ∆s, x i+r+1 ) -u (n∆t, x i+r+1 ) ∆s (21) +v u (n∆t, x i+r+1 ) -u (n∆t, x i+r ) ∆x .

Proposition 2. One has the formula

ε n i = ∆s ∆t ε n i+r+1 .
Proof. The exact solution is constant along the characteristics, that is

u ((n + 1)∆t, x i ) = u (n∆t + ∆s, x i+r+1 ) .
We use this relation in (18) to eliminate u ((n + 1)∆t, x i ). The claim is evident on the resulting quantity.

Proposition 3. Assume that u 0 ∈ W 2,∞ (R). There is a constant C > 0, independent of ∆t, ∆x and v, such that ε n ∞ ≤ Cv u 0 L ∞ (R) ∆x. ( 22
)
Proof. A standard Taylor-Lagrange expansion on the error ε n i gives

ε n i = (∂ t u + v∂ x u) (n∆t, x i ) + ∆s 2 ∂ tt u(n∆t + τ ∆s, x i ) + v∆x 2 ∂ xx u(n∆t, x i + ξ∆x)
and then

ε n ∞ = O ∆s ∂ tt u L ∞ (R) + O v ∆x ∂ xx u L ∞ (R) = O v 2 ∆s ∂ xx u L ∞ (R) + O v ∆x ∂ xx u L ∞ (R) = O (1 + ν)v∆x u 0 L ∞ (R)
thanks to (20). Using [START_REF] Thomée | Stability of difference schemes in the maximum-norm[END_REF] it proves the claim.

An optimal error estimate

The numerical error is e n = (e n i ) i∈Z . with

e n i = u(n∆t, x i ) -u n i . Proposition 4. One has the error estimate for all n such that n∆t ≤ T e n ∞ ≤ CT u 0 L ∞ (R) v∆s∆x ∆t . ( 23
)
The constant C > 0 is independent of ∆t, ∆x and v.

Proof. Let us denote the iteration operator R ν so that ( 13)-( 17) is rewritten as

u n+1 = R ν u n . Thanks to (21) one has u ((n + 1)∆t, •) = R ν u (n∆t, •) + ∆t ε n .
Therefore the error e n satisfies e n+1 = R ν e n + ∆t ε n . Since the numerical scheme satisfies the maximum principle, then R ν ∞ ≤ 1. Therefore

e n+1 ∞ ≤ e n ∞ + ∆t ε n ∞ .
Summation over n yields e n ∞ ≤ n∆t v∆s ∆t ε n ∞ . Using (22) it proves the claim.

Remark 5. This error estimate is now optimal with respect to ∆t. If ∆s = ∆t which means that the characteristics does not go out the first cell, then we recover the first order estimate of convergence characteristics of the upwind scheme. If ∆s << ∆t one can use the bound v∆s ≤ ∆x which is always true and recover the standard error estimate obtained in the numerical analysis of semi-lagrangian numerical methods. We also notice that the error estimate is optimal with respect to v since the error vanishes if v → 0 + .

Numerical scheme

In this section the interpolation function ũn is obtained by a lagrangian interpolation of p + 1 points. Like in previous section, we introduce r ∈ Z such that x i+r ≤ x i -v∆t < x i+r+1 . and we introduce also k ∈ Z: the expression of ũn on [x i+r+1+k-p , x i+r+1+k ] is obtained with a lagrangian interpolation of the values {u i+r+1+k-p , . . . , u i+r+1+k }

ũn (x) = k l=k-p L i+r+l (x)u n i+r+1+l (24)
where the Lagrange polynomials are We define for convenience

L i+r+1+l (x) = k h=k-p,h =l (x-x i+r+1+h ) k h=k-p,h =l (x i+r+l -x i+r+1+h ) . t n (n+1) ∆ x ∆ t t r=-2
α l (ν, k, p) = L i+r+1+l (x i -v∆t) = k h=k-p,h =l (h + ν -1) k h=k-p,h =l (h -l) .
The reduced Courant number is ν = v∆s ∆x ∈]0, 1] and ∆s = (r + 1) ∆x v + ∆t as before. The scheme can be written as

u n+1 i = k l=k-p α l (ν, k, p) u n i+r+1+l , i ∈ Z. ( 25
)
Four parameters characterize this formula: p + 1 is the interpolation order, k determines the local stencil shift and ν is the reduced Courant number; the number of cells that are crossed by the characteristics is (11) is r. The iteration operator is conveniently defined as R ν,k,p so that (25) is equivalent to

u n+1 = R ν,k,p u n . ( 26
)

Stability properties

Introducing the discrete Fourier transform of the vector

u n ūn (ψ) = ∞ j=-∞ e -ijψ u n j , i 2 = -1, (27)
we get from ( 25)

ūn+1 (ψ) = λ ν,k,p (ψ)ū n (ψ), (28)
where

λ ν,k,p (ψ) = k l=k-p α l (ν, k, p) e i (r+1+l)ψ (29)
is the amplification factor for Fourier modes of this scheme. We also introduce the reduced amplification factor is given by

µ ν,k,p (ψ) = k l=k-p α l (ν, k, p) e i lψ , ( 30
)
and we have the relation

λ ν,k,p (ψ) = e i(r+1)ψ µ ν,p,k (ψ). (31)
Proposition 7. The amplification factor satisfies

|λ ν,k,p (ψ)| ≤ 1, ∀ψ ∈ R (32)
if and only if p ∈ {2k, 2k + 1, 2k + 2}.

Many proves are available in the case r + 1 = 0. We refer to the seminal work of Strang in another context [START_REF] Strang | Trigonometric polynomials and difference methods of maximum accuracy[END_REF][START_REF] Iserles | The optimal accuracy of difference schemes[END_REF][START_REF] Iserles | Order stars[END_REF]. See also [START_REF] Després | Uniform asymptotic stability of Strang's explicit compact schemes for linear advection[END_REF][START_REF] Besse | Convergence of classes of high-order semilagrangian schemes for the Vlasov-Poisson system[END_REF][START_REF] Hundsdorfer | Numerical solution of time-dependent advection-diffusion-reaction equations[END_REF]. In consequence one has the following stability property. Proposition 8. Assume that p ∈ {2k, 2k + 1, 2k + 2}. The iteration operator is bounded as an operator over discrete functions with bounded

L 2 norm. More precisely R ν,k,p u 2 ≤ u 2 , ∀u = (u j ) j∈Z . ( 33
)

Interpolation-based error estimate in L 2

We introduce the numerical error of the scheme defined by e n i := u(n∆t, x i )-u n i and e n = (e n i ) i∈Z . One has the standard decomposition

e n+1 = R ν,k,p e n + g n ( 34
)
where g n = (g n i ) is the interpolation error

g n i = u (n∆t, x i -v∆t) - k l=k-p α l (ν, k, p)u (n∆t, x i+r+1+l ) .
We will use the following standard interpolation formula for which we refer the reader to [START_REF] Schumaker | Spline functions: basic theory[END_REF] page 124. Let Q p+1 i be the spline function over the p + 2 points x i -v∆t and

x i+r+1+l for l = k -p ≤ i ≤ k. Q p+1 i
is a piecewise polynomial of degree n and has compact support in ]x i+r+1+k-p , x i+r+1+k [ since x i -v∆t ∈]x i+r+1+k-p , x i+r+1+k [ by definition of r. The interpolation formula writes

g n i = ω i p! ˆxi+r+1+k x i+r+1+k-p Q p+1 i (t)u (p+1) (t)dt (35)
with

ω i = l=i+r+1+k l=i+r+1+k-p (x i -v∆t -x l ).
Other properties of Q p+1 i are ( [START_REF] Schumaker | Spline functions: basic theory[END_REF] page 124) ˆxi+r+1+k

x i+r+1+k-p Q p+1 i (t)dt = 1 p + 1 , ( 36
) and 0 ≤ Q p+1 i (t) ≤ 1 x i+r+1+k -x i+r+1+k-p = 1 p∆x . ( 37
)
Proposition 9. One has the inequality

g n 2 ≤ C k,p ν(1 -ν) u (p+1) 0 L 2 (R) ∆x p+1 ( 38
)
where

C k,p = 1 (p + 1) 1 2 ×          (k+1)!k! (2k+1)! p = 2k + 1, (k)!(k)! (2k)! p = 2k, (k+1)!(k+1)! (2k+2)! p = 2k + 2.
Proof. First we use (19) to eliminate v∆t = v∆s-(r +1)∆x in the interpolation formula (35). We obtain

ω i = k l=k-p (-v∆s -l∆x) = (-1) p+1 ∆x p+1 k l=k-p (l + ν), so |ω i | ≤ ∆x p+1 k l=k-p (l + ν) .
Consider for example the case p = 2k + 1. A rearrangement shows that

k l=k-p (l + ν) = ν(1 -ν) (1 + ν)(2 -ν) ≤1×2 • • • (k + ν)(k + 1 -ν) ≤k×k+1 ≤ ν(1 -ν) (1 × • • • × k) (2 × • • • × (k + 1)) ≤ ν(1 -ν)p!(p + 1) 1 2 C k,p . (39)
A similar trick shows that Π

k l=k-p (l + ν) ≤ ν(1 -ν)p!C k,p (p + 1)
1 2 in all cases. Secondly we set I i =]x i+r+1+k-p , x i+r+1+k [. The Cauchy-Schwarz inequality applied to (35) yields

|g n i | ≤ |ω i | p! Q p+1 i L 2 u (p+1) L 2 (Ii) ≤ |ω i | p! Q p+1 i 1 2 L 1 Q p+1 i 1 2 L ∞ u (p+1) L 2 (Ii) ≤ |ω i | p! 1 p + 1 1 2 1 p∆x 1 2 u (p+1) L 2 (Ii) . Therefore g n 2 = ∆x i |g n i | 2 1 2 ≤ ∆x 1 2 |ω i | p! 1 p + 1 1 2 1 p∆x 1 2 i u (p+1) 2 L 2 (Ii) 1 2 . We notice that i u (p+1) 2 L 2 (Ii) = p u (p+1) 2 L 2 (R) = p u (p+1) 0 2 L 2 (R)
.

Finally we use all these inequalities and obtain the claim (38) after simplification.

Proposition 10. Assume that p ∈ {2k, 2k + 1, 2k + 2}. Then the following error estimate hold for all n such that n∆t ≤ T :

e n 2 ≤ C k,p (1 -ν) v∆s∆x p ∆t T u (p+1) 0 L 2 (R) . ( 40
)
Proof. One gets from (33-34) e n+1 2 ≤ e n 2 + g n 2 , that is, since e 0 = 0,

e n 2 ≤ nC k,p ν(1 -ν) u (p+1) 0 L 2 (R) ∆x p+1 .
It proves the claim after rearrangements.

Remark 11. A crude estimate based on the Stirling formula shows that

C k,p = O 2 -2k
is very small for large p.

Proof. Consider for example the case p = 2k + 1: the constant is

C k,p = 1 2 2k+2 × 2 2k+1
(p + 1)

1 2 2k + 2 k + 1 . ( 41
)
Since the Stirling formula shows that

2 2k+1 2k + 2 k + 1 ≈ Ck 1 2 for some constant C
when k is large, it shows the bound of the remark.

Truncation-based error estimate in L 2

In this section we propose another analysis of the error, which is based on the analysis of the truncation error. The philosophy is rather different. The analysis will fill the gap between interpolation error analysis based on (35) and the Fourier analysis as used in the work of Thomée [START_REF] Thomée | Besov spaces and applications to difference methods for initial value problems[END_REF]. The conclusions will be the same, but with a slightly better constant for large p. This analysis is confirmed with a different approach discussed in the appendix. In both cases the key is the use of advanced formulas which ultimately may be used to provide sharper L 2 estimates for B-splines.

Truncation error

It is convenient to write the scheme in the incremental form

u n+1 j -u n j ∆t + v u n j+ 1 2 -u n j-1 2 ∆x = 0 ( 42
)
where u n j+ 1 2 is the flux. The Fourier symbol of the flux is the function ψ → γ ν,k,p (ψ) such that ūj+ 1 2 (ψ) = γ ν,k,p (ψ)ū j (ψ), where ū is defined by (27). We deduce from the expressions (28) and (42) the following expression for γ ν,k,p : Proposition 12. The Fourier symbol of the flux is

γ ν,k,p (ψ) = λ ν,k,p (ψ) -1 (e iψ -1) v∆t ∆x . ( 43
)
The method we propose here is use this expression of the Fourier symbol of the flux in the truncation error ε n i to obtain the analogous of Proposition 2 for arbitrary order schemes. The truncation error can here be written under the form

ε n i = u((n + 1)∆t, x i ) -u(n∆t, x i ) ∆t + v φ (u(n∆t, •)) (x i ) -φ (u(n∆t, •)) (x i-1 ) ∆x . ( 44
) Let u(•, θ) = ˆR u(•, x) e -iθx dx. ( 45
)
be the Fourier transform u in the x variable of u. Then the Fourier symbol of the operator φ in (44) is given by γ ν,k,p (θ∆x). Moreover, we easily verify, thanks to the relations u(n∆t+∆t, x) = u(n∆t, x-v∆t) and v∆t = v∆s -(r + 1)∆x that u(n∆t + ∆t, θ) = e iθ((r+1)-ν)∆x u(n∆t, θ).

Assuming that u(t, •) ∈ L 2 (R) for all t, we have the following relation

u(•, x) = 1 2π ˆR u(•, θ)e iθx dθ. ( 46
)
Now let express with thanks to the inversion formulae (46) the truncation error :

ε n i = 1 2π ˆR e i((r+1)-ν)θ∆x -1 ∆t + v γ ν,k,p (θ∆x) 1 -e iθ∆x ∆x u(n∆t, θ)e iθxi dθ = 1 2π ˆR e i((r+1)-ν)θ∆x -λ ν,k,p (θ∆x) ∆t u(n∆t, θ)e iθxi dθ = 1 2π ∆s ∆t ˆR e -iνθ∆x -µ ν,k,p (θ∆x) ∆s u(n∆t, θ)e iθxi+r+1 dθ
(47) thanks to (43) and (31). We now recognize at the right side the truncation operator ε n i+r+1 of the semi-lagrangian scheme but with a time step equal to ∆s. At last, we obtain the following Proposition, Proposition 13. One has the formula ε n i = ∆s ∆t ε n i+r+1 for all i and n, and for all coefficients ν, p, k. This expression is the generalization of Proposition 2 at any order.

Bounds

To continue we use a formula which is proved in [START_REF] Després | Uniform asymptotic stability of Strang's explicit compact schemes for linear advection[END_REF]. We consider here the scheme of time step ∆s and the amplification factor µ ν,k,p ; the truncation error corresponding to this scheme is εn i . According to [START_REF] Després | Uniform asymptotic stability of Strang's explicit compact schemes for linear advection[END_REF] (1.5), the amplification factor can be written on the following form :

µ ν,k,p (ψ) = e -iνψ (1 -α(ψ))
where we have set for convenience

α(ψ) = i p+1 α k,p (ν)2 p ˆψ 0 sin p ϕ 2 e i(k-p 2 +ν)(ϕ) dϕ (48)
and

α k,p (ν) = p q=0 (k + ν -q) p! . ( 49
)
Therefore the corresponding truncation error is, thanks to (47)

ε n i = 1 2π 1 ∆s ˆR α(θ∆x) u(n∆t, θ)e iθxi dθ ( 50
)
A direct use of the kernel α is possible using the theory of Thomée [START_REF] Thomée | Stability of difference schemes in the maximum-norm[END_REF][START_REF] Thomée | Besov spaces and applications to difference methods for initial value problems[END_REF]. But it is much more efficient for our purposes to use the following trick. We rewrite first (50) with a backward Fourier transform u(n∆t, θ) = ˆR u(n∆t, x)e -iθx dx.

One obtains

ε n i = 1 2π∆s ˆR u(n∆t, x) ˆR α(θ∆x)e iθ(xi-x) dθ dx. ( 51
)
Then we integrate p + 1 times by parts with respect to x

ε n i = 1 2π∆s ˆR ∂ (p+1) x u(n∆t, x) ˆR α(θ∆x) i p+1 θ p+1 e iθ(xi-x) dθ dx.
A restriction is that we must consider p ≥ 2 so that the internal integral is absolutely convergent. This is not a real restriction since the case p = 1 is elementary and can be treated separately as in Section 3. The kernel in the integral is

β i (x) := ˆR α(θ∆x) i p+1 θ p+1 e iθ(xi-x) dθ. ( 52
)
One has the property :

Lemma 14. The kernel β i has compact support in [x i + (k -p)∆x, x i + k∆x].
Proof. The formula

ε n i = 1 2π∆s ´R ∂ (p+1) x
u(n∆t, x)β i (x)dx shows that the truncation error is a function of ∂ (p+1) x u only, but not of the other derivatives. It is evident that there exists an integral representation

ε n i = ˆxi+k∆x xi+(k-p)∆x K(x)∂ (p+1) x u(n∆t, x)dx in the compact interval [x i + (k -p)∆x, x i + k∆x]
for some kernel K. This is due to two facts: the stencil of the scheme is compact; the order of the scheme is p. Since the function ∂ (p+1) x u is arbitrary in these representation formulas, it implies that

β i (x) = 2π∆s K(x) if x i + (k -p)∆x < x < x i + k∆x, = 0 otherwise.
This is why β i has indeed a compact support

ε n i = 1 2π∆s ˆxi+k∆x xi+(k-p)∆x ∂ (p+1) x u(n∆t, x)β i (x)dx. ( 53
)
See also Remark 17.

Proposition 15. There exists C > 0 such that

ε n 2 ≤ C C k,p p 1 4 (1 -ν)(v∆x p ) u (p+1) 0 L 2 (R) ( 54
)
where C k,p is defined by (39).

Proof. Thanks to (53), we get

|ε n i | ≤ 1 2π∆s β i L 2 (R) ∂ (p+1) x u(n∆t, •) L 2 (xi+(k-p)∆x,xi+k∆x)
.

Since i ∂ (p+1) x u(n∆t, •) 2 L 2 (xi+(k-p)∆x,xi+k∆x) = p ∂ (p+1) x u(n∆t, •) 2 L 2 (R)
one gets that

ε n 2 = ∆x i |ε n i | 2 1 2 ≤ (p∆x) 1 2 2π∆s sup i β i L 2 (R) ∂ (p+1) x u(n∆t, •) L 2 (R)
. Now let us estimate β i L 2 (R) , for i ∈ Z. According to (52) and Proposition 14, we have

β i L 2 (R) = ˆR ˆR α(θ∆x) i p+1 θ p+1 e iθ(xi-x) dθ 2 dx 1/2 ≤ C ˆR |α(θ∆x)| 2 θ 2p+2 dθ 1 2 ≤ C|α k,p (ν)|2 p   ˆR ˆθ∆x 0 | sin p ( ϕ 2 )|dϕ 2 dθ θ 2p+2   1 2 . ( 55
)
thanks to (48). Elementary changes of variable show that

ˆR ˆθ∆x 0 | sin p ( ϕ 2 )|dϕ 2 dθ θ 2p+2 = ∆x 2p+1 2 2p-1 ˆ∞ 0 ˆθ 0 | sin p (ϕ)|dϕ 2 dθ θ 2p+2 .
Using (82) (see Annexe B) we get that

β i L 2 (R) ≤ M ∆x p+ 1 2 |α k,p (ν)| p 1+ 1 4
from which we deduce

ε n 2 ≤ C ∆s |α k,p (ν)| p 1 2 + 1 4 ∆x p+1 u (p+1) 0 L 2 (R)
.

Finally we notice that |α k,p (ν)| = 1 p! k l=k-p
(l + ν) and thanks to estimate (39), we obtain

ε n 2 ≤ C C k,p p 1 4 (1 -ν)ν ∆x p+1 ∆s u (p+1) 0 L 2 (R) . ( 56
)
Since ν = v∆s/∆x, the claim is proved.

Theorem 16. Assume p = 2k, 2k + 1 or 2k + 2. One has the optimal error estimate for all n such that n∆t ≤ T :

e n 2 ≤ D k,p 1 - v∆s ∆x v∆s∆x p ∆t T u (p+1) 0 L 2 (R) (57)
Proof. From Proposition 13 one gets the estimate for the true truncation error

ε n 2 ≤ D k,p (1 -ν) v∆s∆x p ∆t u (p+1) 0 L 2 (R)
.

Since by definition of the truncation error

e n+1 = R ν,k,p e n + ∆tε n (58)
one has the estimate e n+1 2 ≤ e n 2 +∆t ε n 2 . After summation over n time steps (n∆t ≤ T ), one gets the result.

Remark 17. The comparison of ( 34) and ( 58) shows that g n = ∆tε n , that is the interpolation error is proportional to the truncation error. In consequence the kernel Q p+1 i is proportional to the kernel β i . So we understand that the Fourier formula ( 55) is just a more accurate way to bound the L 2 norm of the B-spline, at least more accurate than the crude estimate of

Q p+1 i L ∞
together with

Q p+1 i L 2 ≤ Q p+1 i 1 2 L 1 Q p+1 i 1 2 L ∞ (59)
that was used in the proof of Proposition 9. This is why the constant D k,p is better than the constant C k,p by a factor p -1 4 for large p. More material about the L ∞ norm of Q p+1 i is provided in Appendix A.

L ∞ and L 1 error estimates for odd order schemes

We first have the fundamental result [START_REF] Després | Uniform asymptotic stability of Strang's explicit compact schemes for linear advection[END_REF].

Proposition 18. Assume that p = 2k+1. Then powers of the iteration operator are uniformly bounded in L q norms. More precisely there exists C > 0 such that

(R ν,k,p ) l q ≤ C, ∀l ∈ N, ∀ 1 ≤ q ≤ ∞, ∀v ∆t ∆x ≤ 1. ( 60
)
It must be understood that r = -1.

The shifting technique of semi-lagrangian schemes allows to treat all values of the CFL number, that is all values of r. We obtain Proposition 19. Assume that p = 2k+1. Then powers of the iteration operator are uniformly bounded in L q norms. More precisely there exists C > 0 such that

(R ν,k,p ) l q ≤ C, ∀l ∈ N, ∀ 1 ≤ q ≤ ∞, ∀r. (61)
It is then easy to generalize the approach developed in Proposition 9 to obtain Theorem 20. Assume p = 2k + 1. One has the optimal error estimate for all n such that n∆t ≤ T e n q ≤ E k,p 1 -

v∆s ∆x v∆s∆x p ∆t u (p+1) 0 L q (R) . ( 62
)

Periodic domains

All results of the previous sections hold true in periodic domains. In particular we will use the L 2 error estimate (40) or (57) in the next section.

The Vlasov-Poisson equation

We propose in this section to apply estimates obtained in the previous section to improve the error estimate proposed in [START_REF] Besse | Convergence of classes of high-order semilagrangian schemes for the Vlasov-Poisson system[END_REF] in the context of semi-lagrangian schemes for the Vlasov-Poisson system. The adimensional Vlasov-Poisson equation in one dimension in space and in velocity reads

∂f ∂v + v ∂f ∂x + E(t, x) ∂f ∂v = 0, ( 63
) with ∂E ∂x (t, x) = ˆ+∞ -∞ f (t, x, v)dv -1, (64)
where f (t, x, v) is the distribution function of charged particles (ions or electrons) and E(t, x) the self-consistent electric field. We consider periodic boundary conditions on the variable x ∈ [0, L], that is

f (t, 0, v) = f (t, L, v), v ∈ R, t ≥ 0, (65)
and

E(t, 0) = E(t, L). ( 66
)
In order to have a well-posed problem, we add to equation ( 63)-(66) the zeromean condition on the field

E ˆL 0 E(t, x)dx = 0, t ≥ 0 ( 67
)
and a initial condition

f (0, x, v) = f 0 (x, v), x ∈ L, v ∈ R. (68)
The electric field E is given by

E(t, x) = ˆL 0 K(x, y) ˆR f (t, y, v)dv -1 dy
where K is computed from the kernel of the Poisson equation. We refer to [START_REF] Schaeffer | Global existence of smooth solutions to the Vlasov poisson system in three dimensions Communications in Partial Differential Equations[END_REF], for example, for the existence, uniqueness and regularity of the Vlasov-Poisson system (63)-(68). We will also assume that the initial data has compact support. In consequence the solution has also a compact support in velocity. It is enough to consider the domain v ∈ [-v max , v max ] for t ≤ T . The domain of simulation can be taken as

Ω = [0, L] per × [-v max , v max ] per .

Numerical strategy and convergence

Let N x , N v ∈ N * . We a define numerical approximation f

(n) ∈ R Nx×Nv , by f (n) i,j , i = 0, . . . , N x -1, j = 0, . . . , N v -1, n ∈ N.
The L 2 norm of a discrete function is defined in [START_REF] Després | Uniform asymptotic stability of Strang's explicit compact schemes for linear advection[END_REF]. We introduce the reconstruction operator in the

x direction R x : R Nx×Nv → L ∞ ([0, L]) × R Nv or R x : R Nx → L ∞ ([0, L]) f (•) = R x f.
In practice we use a Lagrangian interpolation of order p + 1 in the x direction as described in (24). We have the same definition for the reconstruction in the v direction R v :

R Nx×Nv → R Nx × L ∞ ([-v max , v max ]) or R v : R Nv → L ∞ ([-v max , v max ]) f (•) = R v f,
with also a Lagrangian reconstruction of order p + 1 but in the v direction. We also define the reconstruction in both directions R x ⊗ R v : R Nx×Nv → L ∞ (Ω m ) which is the "tensorial" product of the one dimensional reconstructions R x and R v : that is we reconstruct first in the v direction and second in the x direction. Notice that the result may be different if one reconstructs first in the x direction and after that in the v direction, unless a compatibility condition is satisfied for commutativity. In the context of this work, such a compatibility condition is not required: the only thing is to guarantee the error estimate (72). For a continuous function g ∈ C(Ω), we define the projected discrete function Πg ∈ R Nx×Nv , by (Πg) i,j = g(x i , v j ), i = 0, . . . , N x -1, j = 0, . . . , N v -1.

Definition 21. We define the discrete transport operator in the x direction based on a semi-lagrangian scheme of order p + 1 by T x,s : R Nx×Nv → R Nx×Nv . Definition 22. We define the discrete transport operator in the v direction based on a semi-lagrangian scheme of order p + 1 by T v,E,s : R Nx×Nv → R Nx×Nv .

A fundamental result deduced from the stability property (33) is

T x,s 2 ≤ 1 and T v,E,s 2 ≤ 1 (69)
for all s and E. This is automatically true is one uses the semi-Lagrangian schemes or Strang's stencil discussed in the first part of this paper. We will also use the exact transport operators T x,s and T v,E,s defined for functions.

Algorithm

The Vlasov-Poisson discrete scheme [START_REF] Falcone | Optimal discretization steps in semi-Lagrangian approximation of first order PDEs[END_REF][START_REF] Pinto | Convergence of an adaptive semi-Lagrangian scheme for the Vlasov-Poisson system[END_REF][START_REF] Besse | Convergence of classes of high-order semilagrangian schemes for the Vlasov-Poisson system[END_REF] reads n) , n ∈ N, f (0) = Πf (0), with the electric field calculated with the following approximation of the exact kernel

f (n+1) = T x,∆t/2 T v,E n h ,∆t T x,∆t/2 f (
E n h (x) = ˆL 0 K(x, y) ˆvmax -vmax R x ⊗ R v T x,∆t/2 f (n) (y, v)dv -1 dy. ( 70
)
Remark 23. It can be checked that ( 70) is equivalent to

E n h (x) = ˆL 0 K(x, y)R x ˆvmax -vmax R v T x,∆t/2 f (n) (•, v)dv -1 dy
which is more convenient for implementation purposes. On can notice that the internal integral is

ρ i = ˆvmax -vmax R v T x,∆t/2 f (n) i,• (v)dv -1 = ∆v Nv-1 j=0 T x,∆t/2 f (n) i,j -1,
which provides an easy way to compute it. Another possibility could be to define the discrete electric field with

E n h (x) = ˆL 0 K(x, y) ˆvmax -vmax T x,∆t/2 R x ⊗ R v f (n) (y, v)dv -1 dy.
This formula seems less evident to implement. For theoretical considerations all formulas are equivalent provided the high order error formula (72) holds.

One notices that E n h may be calculated for any x in the domain, not only at grid points x = x i . We will use this trick in the term ε 4 in the decomposition that follows.

Error decomposition

The numerical error is by definition e

(n) ∈ R Nx×Nv e (n) = Πf (t n ) -f (n) .
We have the following error decomposition

e (n+1) = ε 1 + ε 2 + T x,∆t/2 ε 3 + T x,∆t/2 ε 4 + T x,∆t/2 T v,E n h ,∆t ε 5 , with ε 1 =Πf (t n+1 ) -Π T x,∆t/2 T v,E,∆t T x,∆t/2 f (t n ), ε 2 = Π T x,∆t/2 -T x,∆t/2 Π T v,E,∆t T x,∆t/2 f (t n ). ε 3 =Π T v,E,∆t -T v,E n h ,∆t T x,∆t/2 f (t n ), ε 4 =Π T v,E n h ,∆t T x,∆t/2 f (t n ) -T v,E n h ,∆t Π T x,∆t/2 f (t n ), ε 5 =Π T x,∆t/2 f (t n ) -T x,∆t/2 f (n) ,

Time error of the Strang's splitting

We recall the following result (cf [START_REF] Pinto | Convergence of an adaptive semi-Lagrangian scheme for the Vlasov-Poisson system[END_REF] for example)

Lemma 24. If f (t n ) is bounded in W 1,∞ (Ω), then f (t n+1 ) -T x,∆t/2 T v,E,∆t T x,∆t/2 f (t n ) L ∞ (Ω) ≤ C T ∆t 3 . We then get ε 1 2 ≤ C 1 ∆t 3 , ( 71
)
since the domain is bounded and with periodic conditions, and the "initial data" f (t n ) has compact support inside Ω.

Estimate for the electric field

The third term in the decomposition is

(ε 3 ) i,j =f (t n , x i -(v j -E n (x i )∆t)∆t/2, v j -E n (x i )∆t) -f (t n , x i -(v j -E n h (x i )∆t)∆t/2, v j -E n h (x i )∆t) ,
where we denote E n the electric field calculated from T x,∆t/2 f (t n )

E n (x) := ˆL 0 K(x, y) ˆR T x,∆t/2 f (t n ) (y, v)dv -1 dy
One has

E n (x i ) -E n h (x i ) = ˆL 0 K(x i , y) ˆvmax -vmax g(y, v)dvdy, where g = T x,∆t/2 f (t n ) -R x ⊗ R v Π T x,∆t/2 f (t n ) + R x ⊗ R v (Π T x,∆t/2 -T x,∆t/2 Π)f (t n ) + R x ⊗ R v T x,∆t/2 (Πf (t n ) -f (n) ).
By using ( 13) in [START_REF] Besse | Convergence of classes of high-order semilagrangian schemes for the Vlasov-Poisson system[END_REF], we get the following a priori estimate max i∈{0,...,Nx-1}

|E n (x i ) -E n h (x i )| ≤ C max(∆x, ∆v) p+1 + e (n) 2 . (72) Moreover, since E n ∈ L ∞ (O, L), we get max i∈{0,...,Nx-1} |E n h (x i )| ≤ C (1 + e (n) 2 + max(∆x, ∆v) p+1 ) ≤ C (1 + f (n) 2 + Πf (t n ) 2 ) (73)
and then sup

n≤ T ∆t max i∈{0,...,Nx-1} |E n h (x i )| < +∞ (74)
by using the fact that

Πf (t n ) 2 ≤ 2v max L f (t n ) L ∞ (Ω) ≤ 2v max L f (0) L ∞ (Ω) ,
and, from the stability estimates (69)

f (n) 2 ≤ f (0) 2 ≤ 2v max L f (0) L ∞ (Ω) .

Final proof of Theorem 1

Proof. The triangular inequality implies that

e (n+1) 2 ≤ ε 1 2 + ε 2 2 + ε 3 2 + ε 4 2 + ε 5 2 .
The first term is bounded using (71). The second term is bounded as

ε 2 2 ≤ C 2 ν∆x p+1
using (34)-(38). Using (72) the third term is bounded as

ε 3 2 ≤ C 3 ∆t max(∆x, ∆v) p+1 + e (n) 2 .
The fourth term can also be written

ε 4 = Π T v,E n h ,∆t -T v,E n h ,∆t Π T x,∆t/2 f (t n ) (75)
and then

(ε 4 ) 2 ≤ C 4 τ v ∆v p+1
where τ v is an upper bound of all ν =

|E n h (xi)|∆s ∆v for all possible columns in the domain of computation. One must be careful that ∆s is a function of v, but anyway we have ∆s ≤ ∆t. Therefore, a sharp estimate for τ v is, thanks to (74),

τ v = min sup n≤ T ∆t max i∈{0,...,Nx-1} |E n h (x i )| ∆t ∆v , 1 . (76)
Finally the fifth term can also be written

ε 5 = Π T x,∆t/2 -T x,∆t/2 Π f (t n ) + T x,∆t/2 e (n) , ( 77
)
and then is bounded as

ε 5 2 ≤ Π T x,∆t/2 -T x,∆t/2 Π f (t n ) 2 + e (n) 2
which, using (34)-(38), yields

ε 5 2 ≤ C 5 τ x ∆x p+1 + e (n) 2 .
In this formula τ x is an upper bound of all ν = v∆s ∆x for all possible lines in the domain of computation. Since v∆s ≤ ∆x and ∆s ≤ ∆t, one can take as sharp estimate for τ

x τ x = min v max ∆t ∆x , 1 . ( 78
)
Assuming for simplicity that ∆v = α∆x with a parameter α which is independent of ∆x, one finds that, since ∆t ≤ T, ∆x ≤ L,

∆t max(∆x, ∆v) p+1 ≤ max(α p+1 , 1)∆t∆x p+1 ≤ max(α p+1 , 1)(L + T ) min( ∆t ∆x ,

A More about the L ∞ norm of the B-splines

We give more material about B-splines and explain how to recover an optimal L ∞ norm starting from advanced results discussed in [START_REF] Meinardus | On the Chebyshev Norm of Polynomial B-Splines[END_REF]. It provides an alternative way to analyze (59).

A.1 B-Splines

We consider here the B-spline over the interval [0, 1] as in [START_REF] Meinardus | On the Chebyshev Norm of Polynomial B-Splines[END_REF]. In the even case, we consider the B 2d+1,α spline over the 2d + 2 points

0 < 1 2d < • • • < d 2d ≤ d + α 2d < d + 1 2d < • • • < 2d -1 2d < 1,
and in the odd case, we consider the B 2d+2,α spline over the 2d + 3 points

0 < 1 2d + 1 < • • • < d 2d + 1 ≤ d + α 2d + 1 < d + 1 2d + 1 < • • • < 2d 2d + 1 < 1.
The B-splines are here defined so that ˆ1 0 B 2d+1,α (x)dx = 1 2d + 1 , ˆ1 0 B 2d+2,α (x)dx = 1 2d + 2 .

We then look for sup A. Without restriction, we assume p ≥ 2 so that all integrals are convergent. The last integral is bounded as

B A Lemma

ˆ∞ 0 | sin(ϕ)| p ϕ p dϕ = ˆ1 0 • • • + ˆ∞ 1 • • • ≤ ˆ1 0 | sin(ϕ)| p ϕ p dϕ + 1 p + 1 .
Moreover there exists a constant k > 0 such that sin(ϕ) ≤ ϕ-kϕ 3

Figure 1 :

 1 Figure 1: Stencil of the second order lagrangian (p + 1 = 2) interpolant with a shift of one cell: here r = -2.

Figure 2 :

 2 Figure 2: Example of a semi-lagrangian scheme: here r = -2, k = 1 and p = 3.

0≤α<1B

  2d+1,α L ∞ and sup 0≤α<1 B 2d+2,α L ∞ . We have (seeTheorem 5,[START_REF] Meinardus | On the Chebyshev Norm of Polynomial B-Splines[END_REF])B 2d+1,α (x) = x -d+α 2d 2d B 2d+1,α (x) + B2d (x),

Figure 3 : 2 L 2 L 2 and 3 4π d - 3 / 2

 322332 Figure 3: sup 0≤α<1 B 2d+1,α 2 L , sup 0≤α<1 B 2d+2,α 2 L 2 and

Lemma 26 . 2 dθ θ 2p+2 ≤ 1 p + 1

 26211 There exists a constant M > 0 independent of p ≥ 2 such thatProof. Since | sin(ϕ)| ≤ ϕ for ϕ ≥ 0 one gets ´θ 0 | sin p (ϕ)|dϕ ≤ θ p+1 p+1 . So := ˆ∞ 0 ˆθ 0 | sin p (ϕ)|dϕ ˆ∞ 0 ˆθ 0 | sin(ϕ)| p dϕ dθ θ p+1 . Therefore A ≤ 1 p + 1 ˆ∞ 0 ˆ∞ ϕ 1 θ p+1 dθ | sin(ϕ)| p dϕ = 1 p(p + 1) ˆ∞ 0 | sin(ϕ)| p ϕ p dϕ.

. 1 Theorem 1 .

 11 Consider the Vlasov-Poisson equation

	∂f ∂v	+ v	∂f ∂x	+ E(t, x)	∂f ∂v	= 0,	(1)

  p+1 + ∆t3 1 + e c∆t + e c2∆t + • • • + e c(n-1)∆t .

	and thus
		e (n+1)	2 ≤ e (n)	2 + c min	∆t ∆x	, 1 ∆x p+1 + ∆t e (n)	2 + ∆t 3
			≤ e c∆t e (n)	2 + c min	∆t ∆x	, 1 ∆x p+1 + ∆t 3
	since 1 + c∆t ≤ e c∆t . Therefore after summation
	e (n)	2 ≤e cn∆t e (0)	2
	+ c min , 1 ∆x Since e (0) = 0 and n∆t ≤ T , it yields ∆t ∆x
			e (n)	2 ≤ min	∆t ∆x	, 1 ∆x p+1 + ∆t 3 C ∆t	, C > 0.
							1)∆x p+1 ,

It ends the proof of Theorem 1.

  On the other hand, B p+1,α contains p + 2 points, has it support in [0, 1] and ´1 0 B p+1,α (t)dt = 1 p+2 . We deduce thatp(p + 1)∆xS p+1 (t) = (p + 2)B p+1,α (t),The comparison of inequality (81) with the non optimal formula (37) shows that (81) is better by a factor p -1 2 . Using it in (59), we can improve the estimation of Q p+1 i L 2 by a factor p -1/4 and thus obtain (54) by a different method.

	Remark 25.
	2 Link with Q p+1 i
	Q p+1 i	has it support in [x i+r+1+k-p , x i+r+1+k ], contains p + 2 points and
					ˆxi+r+1+k x i+r+1+k-p	Q p+1 i	(x)dx =	1 p + 1	.
	We then consider S p+1 (t) = Q p+1 i	(x i+r+1+k-p + tp∆x), which has it support in
	[0, 1]. We have			ˆ1 0	S p+1 (t)dt =	1 p(p + 1)∆x	.
	and thus	Q p+1 i	(x) =	p + 2 p(p + 1)∆x	B p+1,α	x -x i+r+1+k-p p∆x	.
	We finally get the optimal bound
	0 ≤ Q p+1 i	(x) ≤	p + 2 p(p + 1)∆x	1 π ˆ∞ -∞	sin t t	p+1	dt ∼	1 ∆x	6 π	p -3/2 .	(81)

Remark 6. The estimate (23) can be enhanced with a term 1 -ν in the right hand side. This is due to the fact that the scheme is exact if ν = 1. This term is visible in the general estimate (40).4 Semi-lagrangian interpolation schemes of order p + 1This section is devoted to the study of arbitrary order semi-lagrangian schemes.

where B2d is the uniform B-spline with 2d + 1 points, and also

where B2d+1 is the uniform B-spline with 2d + 2 points. Now, let x * 2d+1,α such that B 2d+1,α (x * 2d+1,α ) = B 2d+1,α L ∞ . We have B 2d+1,α (x * 2d+1,α ) = 0, and thus B 2d+1,α ∞ = B2d (x * 2d+1,α ). We also have B 2d+1,0 (1/2) = B2d (1/2). We then get

which means that (see proof of Theorem 4 in [START_REF] Meinardus | On the Chebyshev Norm of Polynomial B-Splines[END_REF])

In the even case, we get similarly

As a consequence, we get similar formulas for even and odd cases

and

On Figure 3 L 2 in logarithmic scale increases (in absolute value) from -1.08 to -1.35. We think that the the asymptotic value of the slope is ≈ -1.5. It is reasonable to infer that the asymptotic slope is indeed -1.5 which, if it is true, will prove that (79-80) are optimal.