L Gosse
email: l.gosse@ba.iac.cnr.it

N Vauchelet
email: vauchelet@math.univ-paris13.fr

Some examples of kinetic schemes whose diffusive limit is Il'in's exponential-fitting

Keywords: Discrete-Ordinates (DO), Diffusive approximation, Exponential fitting scheme, Kinetic well-balanced scheme, Run-and-Tumble, Uniformly accurate scheme, Vlasov-Fokker-Planck. 2010 AMS subject classifications: 65M06, 34D15, 76M45, 76R50

This paper is concerned with diffusive approximations of peculiar numerical schemes for several linear (or weakly nonlinear) kinetic models which are motivated by wide-range applications, including radiative transfer or neutron transport, run-and-tumble models of chemotaxis dynamics, and Vlasov-Fokker-Planck plasma modeling. The well-balanced method applied to such kinetic equations leads to time-marching schemes involving a "scattering S-matrix", itself derived from a normal modes decomposition of the stationary solution. One common feature these models share is the type of diffusive approximation: their macroscopic densities solve drift-diffusion systems, for which a distinguished numerical scheme is Il'in/Scharfetter-Gummel's "exponential fitting" discretization. We prove that the well-balanced schemes relax, within a parabolic rescaling, towards the Il'in exponential-fitting discretization by means of an appropriate decomposition of the S-matrix. This is the so-called asymptotic preserving (or uniformly accurate) property.

Introduction and contextualization 1.General setup

Drift-diffusion equations, like (1.3), arise naturally as diffusive approximations of numerous kinetic equations when time t and space x variables are conveniently rescaled. Such parabolic equations, in a context of semiconductor modeling, sparkled the development of so-called "uniformly accurate" (nowadays rephrased "asymptotic-preserving" (AP)) numerical methods: see e.g. [START_REF] Ainsworth | Fundamental systems of numerical schemes for linear convectiondiffusion equations and their relationship to accuracy[END_REF][START_REF] De | A suggested approach to finite difference representation of differential equations[END_REF][START_REF] Gartland | On the uniform convergence of the Scharfetter-Gummel discretization in one dimension[END_REF][START_REF] Ilin | A difference scheme for a differential equation with a small parameter affecting the highest derivative[END_REF][START_REF] Roos | Ten ways to generate the Il'in and related schemes[END_REF][START_REF] Scharfetter | Large signal analysis of a silicon read diode oscillator[END_REF]. A thorough survey of such algorithms is presented in the book [START_REF] Roos | Robust numerical methods for singularly perturbed differential equations. Convection-diffusion-reaction and flow problems; 2nd edition[END_REF]. The "exponentialfitting" Il'in/Scharfetter-Gummel algorithm realizes one of the first well-balanced (WB) schemes for a parabolic equation in divergence form, like drift-diffusion equations, and is uniformly accurate (or AP), too, in the vanishing viscosity limit (contrary to the more standard Crank-Nicolson method). Recently well-balanced numerical methods have been proposed in [24, Part II] to discretize kinetic equations. The formulation of these schemes involves so-called scattering S-matrices built on exponential "Knudsen layers". Accordingly, it is quite natural to wonder how this approach may be adapted to build well-balanced numerical methods which 1. are "uniformly accurate" (or AP) within a diffusive rescaling of variables ? 2. lead asymptotically to an "exponential-fitting" discretization ?

In this paper, we intend to give a positive answer to both these questions for three examples of kinetic equations for which an explicit form of the S-matrix is known. Let f (t, x, v) be a distribution function, depending on time t > 0, position x ∈ R, and velocity v ∈ V : we shall consider,

• a first kinetic model, in parabolic scaling, which reads

ε∂ t f ε + v∂ x f ε = 1 ε 1 -1 T ε (t, x, v ′ )f ε (t, x, v ′ ) dv ′ 2 -T ε (t, x, v)f ε (t, x, v) . (1.1)
When T ε ≡ 1, the well-known conservative radiative transfer equation is recovered, which, as ε → 0, approaches the heat equation. When modeling chemotactic motions of bacteria, equation (1.1) is the so-called Othmer-Alt model [START_REF] Othmer | The diffusion limit of transport equations II: Chemotaxis equations[END_REF]. The tumbling rate T ε describes the response to variations of chemical concentration along a bacteria's path. When the parameter ε → 0, it is now well-established that the macroscopic density ρ :=

• a related model, the Vlasov-Fokker-Planck system, for which the integral collision term is reduced to a diffusion operator. It reads, in parabolic scaling,

ε∂ t f ε + v∂ x f ε + E • ∂ v f ε = 1 ε ∂ v (vf ε + κ∂ v f ε ) . (1.2)
It converges, as ε → 0 towards the drift-diffusion equation, [START_REF] Nieto | High field limit for the Vlasov-Poisson-Fokker-Planck system[END_REF][START_REF] Poupaud | Parabolic limit and stability of the Vlasov-Poisson-Fokker-Planck system[END_REF][START_REF] Wollman | Numerical approximation of the Vlasov-Poisson-Fokker-Planck system in one dimension[END_REF].

Scope and plan of the paper

An object lying at the center of our matters is the so-called "exponential-fit" (Il'in [START_REF] Ilin | A difference scheme for a differential equation with a small parameter affecting the highest derivative[END_REF], Scharfetter-Gummel [START_REF] Scharfetter | Large signal analysis of a silicon read diode oscillator[END_REF], or Chang-Cooper [START_REF] Chang | A Practical Difference Scheme for Fokker-Planck Equations[END_REF]) numerical scheme for 1D drift-diffusion equations, that we briefly recall now. Consider,

∂ t ρ -∂ x (D∂ x ρ -Eρ) = 0. 0 ≤ D, E ∈ R. (1.3)
This equation is discretized by means of a conservative numerical flux:

ρ n+1 j -ρ n j ∆t - F n j+ 1 2 -F n j-1 2 ∆x = 0, (1.4)
where ρ n j is an approximation of ρ(t n , x j ). The flux F n j-1 2 is an approximation of J := D∂ x ρ -Eρ at each interface of the grid x j-1 2

, which is derived by taking advantage of stationary solutions.

J = D∂ x ρ(x) -E ρ(x), ρ(0) = ρ n j-1 , ρ(∆x) = ρ n j .
Then, elementary calculations lead to:

F n j-1 2 := J = E ρ n j-1 -exp(-E ∆x/D) ρ n j 1 -exp(-E ∆x/D) (1.5)
Scheme (1.4)-(1.5) constitutes the classical Il'in/Sharfetter Gummel scheme [START_REF] Gartland | On the uniform convergence of the Scharfetter-Gummel discretization in one dimension[END_REF][START_REF] Ilin | A difference scheme for a differential equation with a small parameter affecting the highest derivative[END_REF][START_REF] Roos | Ten ways to generate the Il'in and related schemes[END_REF][START_REF] Scharfetter | Large signal analysis of a silicon read diode oscillator[END_REF]. Our main result may be formulated as:

Theorem 1.1 Let (2.10)-(2.9) be a numerical scheme relying on a S-matrix for any of the linear 1 + 1 kinetic equations (1.1) and (1.2) in parabolic scaling. Then, for 0 ≤ ε ≪ 1, and uniformly in ∆x > 0,

• each one, among the three considered S-matrices acting in (2.10)-(2.9), admits a decomposition of the type (2.11); such a decomposition yields a Well-Balanced/Asymptotic-Preserving IMEX scheme like (2.12);

• when ε → 0 in the numerical scheme, the corresponding macroscopic densities satisfy the Il'In/Sharfetter-Gummel scheme (1.4)- (1.5).

Obviously, such a general statement contains several former ones, among which the two-stream Goldstein-Taylor model relaxing to the heat equation, [START_REF] Gosse | An asymptotic preserving well-balanced scheme for the hyperbolic heat equation[END_REF], or the so-called "Cattaneo model of chemotaxis", [START_REF] Gosse | Asymptotic-Preserving and Well-Balanced scheme for the 1D Cattaneo model of chemotaxis movement in both hyperbolic and diffusive regimes[END_REF]; some of these results were surveyed in [START_REF] Gosse | Computing Qualitatively Correct Approximations of Balance Laws[END_REF]Part II].Yet, as a guideline for more involved calculations, we first explain in §3 how the limiting process works for a simple two-stream approximation of (1.1), the so-called "Greenberg-Alt" model of chemotaxis [START_REF] Greenberg | Stability results for a diffusion equation with functional shift approximating a chemotaxis model[END_REF]. We mention that for this simple two-velocity model, a numerical scheme formulated in terms of (2.10) with 2 × 2 S-matrices is provided in both [24, page 158] and [START_REF] Gosse | Numerical high-field limits in two-stream kinetic models and 1D aggregation equations[END_REF]Lemma 4.1] (for the purpose of hydrodynamic limits, though). Another type of closely related "diffusive limit" involving a 2 × 2 S-matrix was studied in [START_REF] Gosse | A well-balanced and asymptotic-preserving scheme for the one-dimensional linear Dirac equation[END_REF]. This elementary calculation carried out on the two-stream "Greenberg-Alt" model (3.18) reveals why it is rather natural to expect that a well-balanced algorithm (3.24) (based on stationary solutions) may relax, within a parabolic scaling, toward the exponential-fit scheme (3.21)- (3.22) for the corresponding asymptotic Keller-Segel model. However, as our Theorem 1.1 covers also continuous-velocity models discretized with general quadrature rules, we present in §2 our strategy of proof: in particular, the general scheme involving a scattering matrix is presented in (2.10) and the importance of the decomposition of the scattering matrix (2.11) is emphasized. In §4, such a strategy is applied to the simplest case of continuous equation, namely the "grey radiative transfer" model (4.29). For this system, it is shown in Theorem 4.6 that our numerical scheme relaxes to the finite-difference discretization of the heat equation (4.55). In §5, the case of the Othmer-Alt [START_REF] Othmer | The diffusion limit of transport equations II: Chemotaxis equations[END_REF] model of chemotaxis dynamics (5.56)-(5.57) is handled in a similar manner (at the price of more intricate computations, though), yielding asymptotically the scheme (5.59), this is Theorem 5.5. At last, in §6, the case of a Vlasov-Fokker-Planck model (6.74) is studied, and its asymptotic convergence towards (6.76) is studied.

An essential difference between Vlasov-Fokker-Planck kinetic models and the ones involving an integral collision operator (4.29), (5.56) is that, being exponential polynomials, stationary solutions of (6.74) may not constitute Chebyshev T -systems on v ∈ (0, +∞); definition of Chebyshev T -systems is recalled in §2 below. This drawback has to be compensated by supplementary assumptions on the set of discrete velocities, like (6.84) and (6.85). Accordingly, general properties of eigenfunctions for each of the stationary kinetic models are stated in Appendix, along with a new result on exponential monomials, see [START_REF] Krattenthaler | [END_REF].

Remark 1.2 (Notations) When u ∈ R N and v ∈ R M , the matrix u ⊗ v is an element of M N ×M whose coefficients are (u k v ℓ ).
We will also commonly use the abuse of notations

1 1 + u ⊗ v ∈ M N ×M , with coefficients 1 1 + u k v ℓ k,ℓ .
The present work somehow completes the former ones [START_REF] Gosse | A well-balanced scheme able to cope with hydrodynamic limits for linear kinetic models[END_REF][START_REF] Gosse | Numerical high-field limits in two-stream kinetic models and 1D aggregation equations[END_REF][START_REF] Gosse | Hydrodynamic singular regimes in 1+1 kinetic models and spectral numerical methods[END_REF] where hydrodynamic limits, involving finite-time concentrations, were considered; hereafter, diffusive limits yielding smooth solutions are studied, and general conclusions are identical: the schemes proposed in [24, Part II], involving S-matrices based on stationary solutions, yield more accurate discretizations of the asymptotic regime. Namely, the upwind scheme for hydrodynamic limit (instead of Lax-Friedrichs, [START_REF] James | Numerical methods for one-dimensional aggregation equations[END_REF]), and exponential-fitting in the diffusive one.

2 A general strategy for proving Main Theorem 1.1

= (x 1 , x 2 , ..., x n ) ∈ I n , the family of n vectors (f 1 (X), f 2 (X), ..., f n (X)
) is linearly independent. Equivalently, the determinant never vanishes:

∀ (x 1 , x 2 , . . . , x n ) ∈ I n , x 1 < x 2 < . . . < x n , f 1 (x 1 ) f 2 (x 1 ) • • • f n (x 1 ) f 1 (x 2 ) f 2 (x 2 ) • • • f n (x 2 ) . . . . . . . . . . . . f 1 (x n ) f 2 (x n ) • • • f n (x n ) = 0. (2.6)
Such a family F n constitutes a Chebyshev T -system on the interval I.

The simplest example of T -system on I = R is the monomials family, for which the determinant (2.6) is the well-known Vandermonde determinant.

Definition 2.2 Let F = (f 1 , f 2 , ...) be an infinite sequence of functions, continuous on an interval I ⊂ R: it is said to be a Markov system if, for any n ∈ N, the extracted finite family F n ⊂ F is a Chebyshev T -system.

A standard, yet important result is the following:

Proposition 2.3 Let n ∈ N and F n = (f 1 , f 2 , ..., f n ) be
∀(a 1 , a 2 , • • • , a n ) ∈ R n , I ∋ x → n i=1 a i f i (x), (2.7)
admits at most n -1 real roots on I.

Remark 2.4 Determinants of the type (2.6) are called "alternant determinants", see [START_REF] Vein | Determinants and their applications in Mathematical Physics[END_REF]Chapter 4]. Moreover, the Haar property is closely related to "total positivity" of matrices, see e.g. [START_REF] Gasca | Total Positivity and Its Applications, Series Mathematics and Its Applications[END_REF]. Being the "Hadamard product" the component-wise product of two n × n matrices,

(A • B) 1≤i,j,n := A i,j B i,j ,
Garloff and Wagner, in [21, page 100], explain that the Haar property is not generally preserved by multiplying elements of two T -systems with each other. A first exception is given by two generalized Vandermonde matrices G n = (x α j i ) 1≤i,j,n sharing either the set of points X or the exponents α i 's. A second one is given by non-negative exponential monomials: see Theorem A.1.

General strategy for building AP and WB schemes

Let us first give a Discrete-Ordinates (DO) setup. A spatial domain is gridded with nodes x j = x 0 + j∆x, j ∈ Z and a velocity domain, V , symmetric with respect to 0, by

v k , k = -K, . . . , -1, 1, . . . , K, 0 < v 1 < v 2 < . . . < v K , and v -k = -v k . We denote V = (v 1 , . . . , v K ) ⊤ ∈ R K , V := diag(v 1 , . . . , v K , v 1 , . . . , v K ) ∈ M 2K (R).
Corresponding weights (ω k ) k=1,...,K may be given by a Gauss quadrature, so

V φ(v)dv is approximated by K k=1 ω k (φ(v k ) + φ(-v k )).
The time step will be denoted ∆t > 0.

We consider a kinetic system in parabolic scaling

ε∂ t f + v∂ x f = 1 ε L(f ), t, x, v, ∈ R + * × R × V,
where L is a linear operator, depending on the nature of the given problem. We assume that when ε → 0 the macroscopic density, defined by ρ := V f (v) dv, converges to a solution of a drift-diffusion kind of equation (1.3). Let a discretization of the distribution function at time t n , be (f n j (±v k )) j,k . Our general strategy for building AP and WB schemes follows these steps.

• 1st step. Determination of the S-matrix. In order to build a well-balanced scheme, i.e.

which preserves equilibria, it is important to be able to compute stationary solutions. Let us consider the following stationary problem with incoming boundary conditions on (0, ∆x) for each j = 1, . . . , N x ,

εv∂ x f = L j-1 2 (f ), f (0, V) = f j-1 (V), f (∆x, -V) = f j (-V).
(2.8)

In this system L j-1 2 is a discretization of L on (x j-1 , x j ) such that L j-1 2
is a linear operator. Notice that it is enough to solve the problem with ε = 1 thanks to the change of variable x → x/ε. The unknwon of the problem is the function f and we want to determine the outgoing flux f (∆x, V) f (0, -V) . Since (2.8) is linear, the computation of the outgoing flux involves a so-called scattering matrix S ε j-1/2 defined by

f j-1 2 (V) f j-1 2 (-V) := f (∆x, V) f (0, -V) = S ε j-1 2 f j-1 (V) f j (-V) .
(2.9)

Several S-matrices for various kinetic models are provided in [24, Part II], including the approach based on Case's elementary solutions [START_REF] Aamodt | Useful identities for half-space problems in linear transport theory[END_REF].

• 2nd step. Well-balanced scheme. Once the scattering matrix is known, one may define the well-balanced scheme as (see [START_REF] Gosse | Redheffer products & numerical approximation of currents in one-dimensional semiconductor kinetic models[END_REF])

f n+1 j (V) f n+1 j-1 (-V) = f n j (V) f n j-1 (-V) - ∆t ε∆x V f j (V) -f j-1 2 (V) f j-1 (-V) -f j-1 2 (-V) (2.10)
It verifies the well-balanced property, i.e. stationary states are preserved.

• 3rd step. Asymptotic preserving scheme. Obviously to have an uniformly accurate scheme and to be able to pass to the limit as ε → 0 into the scheme (2.10), we need to treat implicitly terms in 1 ε . To do so, a crucial step is the decomposition

S ε j-1 2 = S 0 j-1 2 + εS 1,ε j-1 2 .
(2.11)

Finally, the scheme (2.10)-(2.9) becomes

f n+1 j (V) f n+1 j-1 (-V) + ∆t ε∆x V f n+1 j (V) f n+1 j-1 (-V) -S 0 j-1 2 f n+1 j-1 (V) f n+1 j (-V) (2.12) = f n j (V) f n j-1 (-V) + ∆t ∆x VS 1,ε j-1 2 f n j-1 (V) f n j (-V)
.

Recall that an approximation of macroscopic density is recovered thanks to the quadrature by

∀j, n, ρ n j := K k=1 ω k (f n j (v k ) + f n j (-v k )).
(2.13)

Strategy for the proof of the main result

The aim of this paper is to prove that, at least for the kinetic systems (1.1) and (1.2) for which scattering matrices are well established, the limit as ε → 0 of scheme (2.12) leads to the Il'In/Sharfetter Gummel scheme (1.4)-(1.5) for the macroscopic density defined in (2.13). A first ingredient in the proof will be to establish that actually the leading order term in the decomposition (2.11) writes with K × K block matrices

S 0 j-1 2 = 0 K S 0 1,j-1 2 S 0 1,j-1 2 0 K and S 1,ε j-1 2 =   S 1,ε 1,j-1 2 S 1,ε 2,j-1 2 S ,ε 3,j-1 2 S 1,ε 4,j-1 2   .
As a consequence, scheme (2.12) can be recast in IMEX (IMplicit-EXplicit [START_REF] Boscarino | Implicit-Explicit Runge-Kutta schemes for hyperbolic systems and kinetic equations in the diffusion limit[END_REF][START_REF] Pareschi | Implicit-Explicit Runge-Kutta schemes for stiff systems of differential equations, Recent Trends in Numerical Analysis[END_REF]) form,

R ε f n+1 j (V) f n+1 j (-V) = εf n j (V) εf n j (-V) + ε∆t ∆x V   S 1,ε 1,j-1 2 f n j-1 (V) + S 1,ε 2,j-1 2 f n j (-V) S 1,ε 3,j+ 1 2 f n j (V) + S 1,ε 4,j+ 1 2 f n j+1 (-V)   , (2.14)
for a matrix R ε (which details depend on each kinetic model) given by,

R ε = εI R 2K + ∆t ∆x V I K -S 0 1,j-1 2 -S 0 1,j+ 1 2 I K -→ ε→0 R 0 := ∆t ∆x V I K -S 0 1,j-1 2 -S 0 1,j+ 1 2 I K .
To ensure global solvability of such a scheme in the diffusive limit, we shall argue according to the analysis performed for continuous equations, see for instance [START_REF] Jüngel | Transport Equations for Semiconductors[END_REF]Chap. 5], namely by proving that:

• Ker(R 0 ) is a vectorial line, closely related to Maxwellian distributions;

• its range is an hyperplane which contains all kinetic distributions with vanishing macroscopic densities (i.e. null moments of order zero).

These conditions are discrete analogues of the Fredholm alternative holding for continuous limits. Yet, we write that kinetic densities have a Hilbert expansion, f = f 0 + εf 1 + . . .. Injecting into (2.14), we first deduce by identifying the term at order 0 in ε that f 0 ∈ Ker(R 0 ). Identifying the terms at order 1 in ε, we get

R 0 {f 1 } n+1 j (V) {f 1 } n+1 j (V) = ({f 0 } n j -{f 0 } n+1 j )(V) ({f 0 } n j -{f 0 } n+1 j )(-V) (2.15) + ∆t ∆x V   S 1,0 1,j-1 2 {f 0 } n j-1 (V) + S 1,0 2,j-1 2 {f 0 } n j (-V) S 1,0 3,j+ 1 2 {f 0 } n j (V) + S 1,0 4,j+ 1 2 {f 0 } n j+1 (-V)   .
This equation admits a solution iff the right hand side belongs to the range of R 0 . Taking moments of order zero allows us to deduce the asymptotic discretization governing the numerical macroscopic density.

To summarize, the main issues are:

• to compute the scattering matrix S and determine its decomposition (2.11);

• to study the kernel and the range of R 0 ;

• to establish that taking the moment of order zero of (2.15) leads to Il'In/Sharfetter Gummel scheme for macroscopic density.

Finally, for future use, we recall the following result:

Lemma 2.5 (Lemma 3.1 in [START_REF] Calvez | Traveling chemotactic aggregates at mesoscopic scale and bi-stability[END_REF], Prop. 1 in [START_REF] Gosse | Redheffer products & numerical approximation of currents in one-dimensional semiconductor kinetic models[END_REF]) Let Γ be the following diagonal matrix,

Γ = diag(ω k v k ) k=1,...,K 0 K 0 K diag(ω k v k ) k=1,...,K ,
then, being given a nonnegative initial data f 0 j (±V) j , the scheme (2.10) preserves both non-negativity and the (discrete) L 1 norm of f n j (±V) as soon as

max(V)∆t ≤ ε∆x, CFL condition, (2.16)
∀ j, ΓS n j-1/2 Γ -1 is left-stochastic (each column summing to 1).

(2.17)

Moreover, it preserves its L ∞ norm as well if S n j-1/2 is right-stochastic (each row summing to 1).

The two-stream Greenberg-Alt kinetic model

For the sake of simplicity, we first start our exposition by a very simple model consisting in a two-velocity kinetic model.

Diffusive limit of the continuous system

The simplest two-velocity kinetic model describing the motion of bacteria by chemotaxis was given in [START_REF] Greenberg | Stability results for a diffusion equation with functional shift approximating a chemotaxis model[END_REF]. Let f + (f -) be the distribution function of right-moving (left-moving) bacteria, the following system describes their motion governed by a run and tumble process (see also e.g. [START_REF] Chalub | Kinetic models for chemotaxis and their drift-diffusion limits[END_REF][START_REF] Dolak | Kinetic models for chemotaxis: Hydrodynamic limits and spatiotemporal mechanisms[END_REF])

ε∂ t f ± ± ∂ x f ± = ± 1 2ε 1 + εφ(∂ x S))f --(1 -εφ(∂ x S))f + . (3.18)
The quantity S is the chemoattractant concentration, and solves

-∂ xx S + S = ρ. (3.19)
Macroscopic quantities being ρ = f + + f -(density) and J = 1 ε (f +f -) (current), adding and subtracting former equations yields,

∂ t ρ + ∂ x J = 0, ε 2 ∂ t J + ∂ x ρ = φ(∂ x S)ρ -J.
Letting formally ε → 0 in the second equation,

J = φ(∂ x S)ρ -∂ x ρ,
and yields the well-known Keller-Segel, system when coupled with (3.19),

∂ t ρ -∂ xx ρ + ∂ x φ(∂ x S)ρ = 0. (3.20)

Exponential-fit scheme for Keller-Segel equation

We first recall the Il'in/Sharfetter Gummel scheme for system (3.20), which writes under the form (1.3) with D = 1 and E = φ(∂ x S). Assuming that approximations (ρ n j ) j of ρ(t n , x j ) and (S n j ) j of S(t n , x j ) are available, we denote ∀j,

φ n j-1/2 = φ S n j -S n j-1 ∆x .
Then, the Il'In/Sharfetter Gummel scheme (1.4)-(1.5) reads, in this framework,

J n j-1/2 = φ n j-1/2 1 -e -φ n j-1/2 ∆x ρ n j-1 -e -φ n j-1/2 ∆x ρ n j , (3.21)
which gives (constant) currents. Then densities are updated by,

ρ n+1 j = ρ n j + ∆t ∆x J n j-1/2 -J n j+1/2 . (3.22)

Asymptotic preserving and well-balanced scheme

We follow the strategy proposed in §2.2. Let us assume that (f n,+ j , f n,- j ) j are known at time t n along with an approximation S n j of S(t n , x j ), giving φ n j-1/2 = φ(

S n j -S n j-1 ∆x
).

• 1st step. Scattering matrix. It is computed by solving the stationary system in (x j-1 , x j ) with incoming boundary conditions,

   ∂ x f ± = 1 2ε (1 + εφ n j-1/2 )f --(1 -εφ n j-1/2 )f + f + (x j-1 ) = f + j-1 ; f -(x j ) = f - j .
(3.23)

The unknown for this problem are interface values f

+ j-1/2 := f + (x j ) and f - j-1/2 := f -(x j-1 ).
System (3.23) may be solved exactly by, first, subtracting both equations,

∂ x (f + -f -) = 0, J := 1 ε (f + -f -) is constant in (x j-1 , x j ).
and then, adding them, so that by denoting ρ = f

+ + f -, ∂ x ρ = φ n j-1/2 ρ -J, e -φ n j-1/2 ∆x ρ j -ρ j-1 = J e -φ∆x -1 φ .
Thus the system to be solved is

f + j-1 -f - j-1/2 = f + j-1/2 -f - j e -φ n j-1/2 ∆x (f - j + f + j-1/2 ) -(f + j-1 + f - j-1/2 ) = e -φ n j-1/2 ∆x -1 εφ n j-1/2 (f + j-1 -f - j-1/2 ).
After easy computations,

f - j-1/2 =f + j-1 + 2εφ n j-1/2 f + j-1 -e -φ n j-1/2 ∆x f - j e -φ n j-1/2 ∆x -1 -εφ n j-1/2 (1 + e -φ n j-1/2 ∆x ) , f + j-1/2 =f - j - 2εφ n j-1/2 f + j-1 -e -φ n j-1/2 ∆x f - j e -φ n j-1/2 ∆x -1 -εφ n j-1/2 (1 + e -φ n j-1/2 ∆x )
.

Following [24, pp.157-58] (or [START_REF] Gosse | A well-balanced scheme able to cope with hydrodynamic limits for linear kinetic models[END_REF]), one may rewrites the latter system by defining a S-matrix,

f + j-1/2 f - j-1/2 = S n j-1/2 f + j-1 f - j , which (with shorthand notation E = e -φ n j-1/2 ∆x ) reads S n j-1/2 =    -2εφ n j-1/2 E-1-εφ n j-1/2 (1+E) 1 + 2εφ n j-1/2 E E-1-εφ n j-1/2 (1+E) 1 + 2εφ n j-1/2 E-1-εφ n j-1/2 (1+E) -2εφ n j-1/2 E E-1-εφ n j-1/2 (1+E)    .
Notice that it is clearly left-stochastic.

• 2nd step. WB scheme. Following [START_REF] Gosse | An asymptotic preserving well-balanced scheme for the hyperbolic heat equation[END_REF], we consider

f +,n+1 j = f +,n j -∆t ε∆x (f +,n+1 j -f + j-1/2 ) f -,n+1 j-1 = f -,n j-1 -∆t ε∆x (f -,n+1 j-1 -f - j-1/2 ).
(3.24)

• 3rd step. AP scheme. For such a simple case, the decomposition (2.11) is straightforward:

S n j-1/2 = 0 1 1 0 + 2ε • φ n j-1/2 E -1 -εφ n j-1/2 (1 + E) -1 E 1 -E . (3.25)
By treating implicitly the first (stiff) term and plugging into the scheme (2.10),

             f +,n+1 j = f +,n j -∆t ε∆x (f +,n+1 j -f -,n+1 j ) + ∆t ∆x J n j-1/2 f -,n+1 j = f -,n j -∆t ε∆x (f -,n+1 j -f +,n+1 j ) -∆t ∆x J n j+1/2 J n j-1/2 = -2φ n j-1/2 f n,+ j-1 -e -φ n j-1/2 ∆x f n,- j e -φ n j-1/2 ∆x -1-εφ n j-1/2 (1+e -φ n j-1/2 ∆x ) . (3.26)
For this simple case, the limit ε → 0 may be performed easily. Adding the first two equations in (3.26), we obtain,

ρ n j = f +,n j + f -,n j , ρ n+1 j = ρ n j - ∆t ∆x (J n j+1/2 -J n j-1/2 ). (3.27)
Letting then ε → 0, we deduce easily from (3.26) that

J n j-1/2 → -2φ n j-1/2 e -φ n j-1/2 ∆x -1 f n,+ j-1 -e -φ n j-1/2 ∆x f n,- j .
Yet, multiplying the first equations of (3.26) by ε and letting ε → 0, at the limit the relation

f +,n+1 j = f -,n+1
j is enforced. Accordingly, the current rewrites

J n j-1/2 → -φ n j-1/2 e -φ n j-1/2 ∆x -1 ρ n j-1 -e -φ n j-1/2 ∆x ρ n j . (3.28)
Injecting (3.28) into (3.27), Il'in scheme for Keller-Segel system (3.20) is found.

Remark 3.1 All computations are explicit, so there is no real need of R ε when passing to the limit. From (3.26), we get the matrix already met in [29, §3],

R ε = ε + ∆t ∆x -∆t ∆x -∆t ∆x ε + ∆t ∆x -→ ε→0 R 0 = ∆t ∆x 1 -1 -1 1 . Thus, Ker(R 0 ) =Vect{(1, 1) ⊤ } and (u, v) ⊤ ∈ Im(R 0 ) iff u + v = 0.
4 Heat equation as a diffusive limit of radiative transfer

Let us consider the kinetic model (1.1) in the simple case T ε ≡ 1:

ε∂ t f + v∂ x f = 1 ε 1 -1 f (t, x, v ′ ) dv ′ 2 -f . (4.29)
The macroscopic density, at the limit ε → 0, satisfies the heat equation

∂ t ρ - 1 3 ∂ xx ρ = 0, ρ(t, x) = 1 -1 f (t, x, v)dv. (4.30)
Equation (4.29) is usually referred to as to "gray radiative transfer"; more specific models can be drawn by replacing the uniform integral kernel 1 2 by an even, nonnegative, function of the velocity variable, 0

≤ K(v) = K(-v), ε∂ t f + v∂ x f = 1 ε K(v) 1 -1 f (t, x, v ′ )dv ′ -f , 1 -1 K(v) = 1.
For this system, the velocity domain being V = (-1, 1), we assume that the set (ω k , v k ) introduced in §2.2 satisfies the mild restrictions,

K k=1 ω k = 1, K k=1 ω k v 2 k = 1 3 , v -k = -v k . (4.31)

Solving the stationary problem

Let us first investigate some useful properties on the stationary equation for general nonnegative kernel T :

Proposition 4.1 Let f (x, v) = exp(-λx)φ λ (v) be a separated-variables solution of the stationary equation v ∂ x f = 1 -1 T (v ′ ) f (x, v ′ )dν(v ′ ) -T (v) f , v ∈ (-1, 1), (4.32)
where ν is a probability measure on (-1, 1). Then, we have

φ λ (v) = 1 T (v) -λv , 1 -1 T (v) T (v) -λv dν(v) = 1. (4.33)
Moreover, the following orthogonality relation holds

∀λ = µ, 1 -1 v φ λ (v)φ µ (v) T (v) dν(v) = 0. (4.34)
In particular, for λ = 0,

φ 0 (v) = 1 T (v) , ∀µ = 0, J µ = 1 -1 v φ µ (v) dν(v) = 0.
The λ are usually called the "eigenvalues" and the corresponding φ λ are the "Case's eigenfunctions".

∀ λ = µ, N k=-N ω k v k φ λ (v k )φ µ (v k )T (v k ) = 0. (4.35)
In particular, for λ = 0,

∀ µ = 0, N k=1 ω k v k (φ µ (v k ) -φ µ (-v k )) = 0. (4.36) Proof. Inserting the ansatz f (x, v) = exp(-λx)φ λ (v) into (4.32) implies (T (v) -λv)φ λ (v) = 1 -1 T (v)φ λ (v) dν(v).
By linearity, φ λ is defined up to a constant, we may fix this constant by imposing

1 -1 T (v)φ λ (v) dν(v) = 1.
It gives the relations (4.33). Then, for two eigenvalues λ = µ, we have

(T (v) -λv)φ λ (v) = 1; (T (v) -µv)φ µ (v) = 1.
We multiply the first identity by T (v)φ µ (v), the second by T (v)φ λ (v), and integrate over dν(v), we obtain after subtracting the resulting identities

(λ -µ) 1 -1 vφ λ (v)φ µ (v)T (v) dν(v) = 0.
We deduce the orthogonality relation in (4.34).

The scattering matrix and its decomposition (2.11)

Following the strategy proposed in Section 2.2, we first determine the scattering matrix and its expansion (2.11). The stationary problem with incoming boundary data reads, for an index j,

εv∂ x f = 1 2 1 -1 f (v ′ )dv ′ -f , on (0, ∆x), (4.37)
f (0, v) = f j-1 (v), f (∆x, -v) = f j (-v), (4.38)
where f j-1 (|v|), f j (-|v|) are "incoming values", and outgoing ones read:

f j-1/2 (|v|) = f (∆x, |v|), f j-1/2 (-|v|) = f (0, -|v|).
In view of Proposition 4.1, equation (4.33) in the particular case T = 1, the discrete eigenelements are given by 1 2

K k=1 ω k 1 1 -λv k + 1 1 + λv k = 1, φ λ (v) = 1 1 -λv .
Clearly, we have that if λ is an eigenvalue, then -λ is also an eigenvalue and φ -λ (v) = φ λ (-v). We observe also that the eigenvalue λ = 0 is double; the eigenvectors are 1 and xεv. (Indeed we verify easily that such functions solve (4.37)). Therefore, a quite general stationary solution is obtained by truncating to the first 2K eigenmodes, (see e.g. [START_REF] Aamodt | Useful identities for half-space problems in linear transport theory[END_REF][START_REF] Cercignani | Mathematical methods in kinetic theory[END_REF][START_REF] Gosse | Computing Qualitatively Correct Approximations of Balance Laws[END_REF])

f (x, v) = a 0 + b 0 (x -εv) + K-1 ℓ=1 a ℓ e -λ ℓ x/ε 1 -λ ℓ v + b ℓ e λ ℓ (x-∆x)/ε 1 + λ ℓ v , λ ℓ ≥ 0.
We denote the vector of so-called "normal modes", 0 ≤ λ := λ 1 , . . . , λ K-1 ⊤ , and the matrix of "Case's eigenfunction"

Φ ± (x) = e -λ x/ε 1 ∓ V ⊗ λ 1 R K e λ (x-∆x)/ε 1 ± V ⊗ λ x1 R K ∓ εV ∈ M K×2K (R), such that f (x, ±V) = Φ ± (x) a b , with the notations a = (a 1 , . . . , a K-1 , a 0 ) ⊤ and b = (b 1 , . . . , b K-1 , b 0 ) ⊤ .
(We recall that the notation

1 1+V⊗λ denotes the matrix in M K (R) whose coefficients are ( 1 1+v k λ ℓ ) k,ℓ
). Thanks to the boundary conditions (4.38) we have

f (0, V) f (∆x, -V) = Φ + (0) Φ -(∆x) a b , f (∆x, V) f (0, -V) = Φ + (∆x) Φ -(0) a b .
Thus, we deduce (see also [24, Chap. 9 pp.175-176]), that the solution of (4.37)-(4.38) is expressed as,

f j-1/2 (V) f j-1/2 (-V) = M ε M -1 ε f j-1 (V) f j (-V) ,
where

M ε := Φ + (∆x) Φ -(0) = e -λ ∆x/ε 1-V⊗λ 1 R K 1 1+V⊗λ ∆x1 R K -εV 1 1+V⊗λ 1 R K e -λ ∆x/ε 1-V⊗λ εV , (4.39)
M ε := Φ + (0) Φ -(∆x) = 1 1-V⊗λ 1 R K e -λ ∆x/ε 1+V⊗λ -εV e -λ ∆x/ε 1+V⊗λ 1 R K 1 1-V⊗λ ∆x1 R K + εV . (4.40)
Thus we have obtained the first statement of the following Proposition:

Proposition 4.3
The scattering matrix for the radiative transfer system is given by

S ε = M ε M -1 ε , (independent of j) (4.41)
where M ε and M ε are given in (4.39)-(4.40). It admits the decomposition

S ε = 0 K I K -ζγ I K -ζγ 0 K + εB ε , where B ε := 1 ε (A ε M -1 ε -A 0 M -1 0 ), (4.42)
with

A ε = M ε - 0 K I K I K 0 K M ε , ζ = 1 1 -V ⊗ λ - 1 1 + V ⊗ λ ∈ M K×K-1 (R),
and

γ ∈ M K-1×K (R), β ⊤ ∈ R K are such that γ 1 1 -V ⊗ λ = I K-1 , γ1 R K = 0 R K-1 , (4.43)
β ⊤ 1 1 -V ⊗ λ = 0 ⊤ R K-1 , β ⊤ 1 R K = 1. (4.44) Remark 4.4
The existence of γ and β is provided by Proposition A.1 (iii) in Appendix. It's good to have an "intuitive idea" of the nature of ζ, γ and β:

• First, ζ : R K-1 → R K converts a set of K -1 spectral coefficients into the restriction to v ∈ V of a kinetic density F (v) having a specific property, ∀v k ∈ V, F (v k ) = -F (-v k ), so K k=1 ω k (F (v k ) + F (-v k )) = 0.
• Then, γ : R K → R K-1 recovers K -1 spectral coefficients corresponding to "damped modes" (Knudsen layers) out of K samples of any kinetic density F (v k ). More precisely, if we consider a given steady kinetic density,

∀v k ∈ V, G(v k ) = a 0 + K-1 ℓ=1 a ℓ 1 -v k λ ℓ , then γ [G(V)] = (a 1 , a 2 , ..., a K-1 ): this is meaningful for ε ≪ 1.
• Oppositely, β detects the "Maxwellian part" in the decomposition of G, that is, the a 0 coefficient, so that β ⊤ [G(V)] = a 0 .

Yet, consider the product ζγ : R K → R K , applied to G(V). It produces,

(ζγ[G])(v k ) = K-1 ℓ=1 a ℓ 1 1 -v k λ ℓ - 1 1 + v k λ ℓ , k = 1, ..., K, so that, ∀v k ∈ V, I K -(ζγ) [G](v k ) = a 0 + K-1 ℓ=1 a ℓ 1 + v k λ ℓ .
Proof.

1. Based on the simple case (3.25), we define A ε such that,

∀ε > 0, S ε = 0 K I K I K 0 K + A ε M -1 ε ,
which clearly yields:

A ε = M ε - 0 K I K I K 0 K M ε = ζe -λ ∆x/ε 0 R K -ζ -2Vε -ζ 0 R K ζe -λ ∆x/ε 2Vε , As ε → 0, A ε → A 0 = 0 K (-ζ 0 R K ) (-ζ 0 R K ) 0 K , (4.45)
along with,

M ε → M 0 = 1 1-V⊗λ 1 R K 0 K×(K-1) 0 K 0 K×(K-1) 1 R K 1 1-V⊗λ ∆x1 R K := M 01 0 K (0 K×(K-1) 1 R K ) M 02 ,
with the notation

M 01 = 1 1 -V ⊗ λ 1 R K , M 02 = 1 1 -V ⊗ λ ∆x 1 R K .
2. We now intend to compute the following inverse,

M -1 0 = M -1 01 0 K -M -1 02 (0 K×(K-1) 1 R K )M -1 01 M -1 02 , (4.46)
where

M -1 01 = γ β ⊤ ∈ M K,K (R), M -1 02 = γ 1 ∆x β ⊤ ∈ M K,K (R)
, being γ a matrix of size (K -1) × K and β an element in R K , such that (4.43) and (4.44) hold. From the expression (4.45) and (4.46), it comes:

A 0 M -1 0 = (-ζ 0 R K )M -1 02 (0 K×(K-1) 1 R K )M -1 01 (-ζ 0 R K )M -1 02 (-ζ 0 R K )M -1 01 0 K • We first observe that (-ζ 0 R K ) γ β ⊤ = (-ζ 0 R K ) γ 1 ∆x β ⊤ = -ζγ.
• Then, using the second identity in (4.43),

A 0 M -1 0 = 0 K -ζγ -ζγ 0 K .
Finally, we reach the decomposition (4.42).

Lemma 4.5 With identical notation as Proposition 4.3, as ε → 0,

B ε = 1 ∆x (2I K -ζγ)Vβ ⊤ -(2I K -ζγ)Vβ ⊤ -(2I K -ζγ)Vβ ⊤ (2I K -ζγ)Vβ ⊤ + o(1). (4.47)
Proof. From (4.42), and using the form,

M ε = M 0 + 0 K×(K-1) 0 R K e -λ ∆x/ε 1+V⊗λ -εV e -λ ∆x/ε 1+V⊗λ 0 R K 0 K×(K-1) εV ,
we get

B ε = 1 ε (A ε -A 0 M -1 0 M ε )M -1 ε = 1 ε (A ε -A 0 )M -1 ε + ζγ( 1 1+V⊗λ )δ ε 0 R K 0 K×(K-1) ζγV 0 K ζγ( 1 1+V⊗λ )δ ε -ζγV M -1 ε ,
where δ ε = 1 ε e -λ ∆x/ε → 0 as ε → 0. Then,

B ε = ζδ ε + ζγ( 1 1+V⊗λ )δ ε 0 R K 0 K×(K-1) (-2I K + ζγ)V 0 K×(K-1) 0 R K ζδ ε + ζγ( 1 1+V⊗λ )δ ε (2I K -ζγ)V M -1 ε . (4.48)
As ε → 0, we get from (4.48)

B ε = 1 ε (A ε M -1 ε -A 0 M -1 0 ) -→ ε→0 B 0 := 0 K×(2K-1) (-2I K + ζγ)V 0 K×(2K-1) (2I K -ζγ)V M -1 0 .
With the expression of the inverse of M 0 in (4.46), we get

B 0 = 1 ∆x (2I K -ζγ)Vβ ⊤ 1 R K β ⊤ -(2I K -ζγ)Vβ ⊤ -(2I K -ζγ)Vβ ⊤ 1 R K β ⊤ (2I K -ζγ)Vβ ⊤ .
Thanks to (4.44), we have β ⊤ 1 R K = 1 and we are done.

Emergence of the macroscopic discretization

Thanks to Proposition 4.3, we have at hand the expansion of the scattering matrix with respect to ε; then, we may write the final scheme (2.12) for the radiative transfer equation,

f n+1 j (V) f n+1 j-1 (-V) + ∆t ε∆x V f n+1 j (V) -(I K -ζγ)f n+1 j (-V) f n+1 j-1 (-V) -(I K -ζγ)f n+1 j-1 (V) = f n j (V) f n j-1 (-V) + ∆t ∆x VB ε f n j-1 (V) f n j (-V) . ( 4
1 ε R ε f n+1 j (V) f n+1 j (-V) = f n j (V) f n j (-V)) + ∆t ∆x V B 1ε f n j-1 (V) + B 2ε f n j (-V) B 3ε f n j (V) + B 4ε f n j+1 (-V) , (4.50)
where we use the notation

B ε = B 1ε B 2ε B 3ε B 4ε
, and

R ε = εI 2N + ∆t ∆x V I K ζγ -I K ζγ -I K I K . (4.51)
Solving this system amounts to the inversion of matrix R ε . By construction, this scheme satisfies the well-balanced property and is asymptotic preserving.

Proposition 4.7 Assume that for all discrete eigenvalues λ i , i = 1, . . . , K -1, identity (4.36) holds, i.e.

K k=1 ω k v k 1 1 -v k λ i - 1 1 + v k λ i = 0.
Then the scheme (4.50) is mass-preserving, uniformly when ε ≪ 1.

Proof. It is mostly a consequence of Lemma 2.5 and Proposition 4.3. Indeed, by Lemma 2.5, ΓS ε Γ -1 needs to be left stochastic for any ε > 0, which implies (by Proposition 4.3) that, for ε ≪ 1, the matrix

Γ 0 K I K -ζγ I K -ζγ 0 K Γ -1 is left stochastic.
Hence, for all j = 1, . . . , K, the sum of the j th column is equal to 1, that is K k=1 ω k v k (ζγ) kj = 0. Accordingly, it comes (by general assumptions) that,

K k=1 ω k v k (ζγ) kj = K k=1 K-1 i=1 ω k v k 1 1 -v k λ i - 1 1 + v k λ i γ ij = 0.

Consistency with the diffusive limit

This subsection is devoted to the proof of Theorem 4.6.

Proof of Theorem 4.6 When ε → 0, we get from (4.48)

B ε -→ B 0 = 1 ∆x (2I K -ζγ)Vβ ⊤ -(2I K -ζγ)Vβ ⊤ -(2I K -ζγ)Vβ ⊤ (2I K -ζγ)Vβ ⊤ .
Moreover, with (4.51),

R ε = R 0 + εI 2K , where R 0 := ∆t ∆x V I K ζγ -I K ζγ -I K I K .
Assuming that f admits a Hilbert expansion f = f 0 + εf 1 + o(ε). We inject into scheme (4.50).

Then by identifying the term in power of ε, we get

R 0 {f 0 } n+1 j (V) {f 0 } n+1 j (-V) = 0, (4.52)
and at order 0 in ε,

R 0 {f 1 } n+1 j (V) {f 1 } n+1 j (-V) = {f 0 } n j (V) -{f 0 } n+1 j (V) {f 0 } n j (-V) -{f 0 } n+1 j (-V) (4.53) + ∆t 2∆x 2 V (2I K -ζγ)Vβ ⊤ ({f 0 } n j-1 (V) -{f 0 } n j (-V)) (2I K -ζγ)Vβ ⊤ ({f 0 } n j+1 (V) -{f 0 } n j (-V))
.

We will make use of Lemma A.2 in the Appendix, which may be applied with µ = λ since (4.36) holds, it gives:

• Ker(R 0 ) = Span(1 R 2K ), • Im(R 0 ) = Z = (Z 1 Z 2 ) ⊤ , Z i ∈ R K such that K k=1 ω k (Z 1k + Z 2k ) = 0 .
Roughly speaking, the range of R 0 is an hyperplane containing kinetic densities which moment of order zero vanishes. Its kernel is the (one-dimensional) vectorial line of constant kinetic densities, which are the Maxwellians for (4.29). We deduce from (4.52)

{f 0 } n+1 j (±V) = ρ n+1 j 2 1 R K , ρ n j = K k=1 ω k (f n j (v k ) + f n j (-v k )). (4.54)
Then, injecting (4.54) into (4.53) and using (4.44), we obtain

R 0 {f 1 } n+1 j (V) {f 1 } n+1 j (-V) = 1 2 (ρ n j -ρ n+1 j )1 R K (ρ n j -ρ n+1 j )1 R K + ∆t ∆x 2 V (2I K -ζγ)V(ρ n j-1 -ρ n j ) (2I K -ζγ)V(ρ n j+1 -ρ n j )
.

Moreover, by definition of R 0 , we have

1 2∆x R 0 V(ρ n j+1 -ρ n j ) V(ρ n j-1 -ρ n j ) = ∆t 2∆x 2 V V(ρ n j+1 -ρ n j ) + (ζγ -I K )V(ρ n j-1 -ρ n j ) (ζγ -I K )V(ρ n j+1 -ρ n j ) + V(ρ n j-1 -ρ n j )
.

Therefore, adding the last two equalities, we deduce that

R 0 {f 1 } n+1 j (V) -1 2∆x V(ρ n j+1 -ρ n j ) {f 1 } n+1 j (-V) -1 2∆x V(ρ n j-1 -ρ n j ) = 1 2 (ρ n j -ρ n+1 j )1 R 2N + ∆t 2∆x 2 V 2 (ρ n j+1 + ρ n j-1 -2ρ n j ).
A solution exists iff the right hand side belongs to Im(R 0 ), so by Lemma A.2,

0 = ρ n j -ρ n+1 j + ∆t ∆x 2 (ρ n j-1 -2ρ n j + ρ n j+1 ) K k=1 ω k v 2 k . (4.55)
We conclude the proof thanks to (4.31).

5 Othmer-Alt model for one-dimensional chemotaxis

Continuous diffusive limit toward Keller-Segel

We consider now a model of chemotaxis in parabolic scaling, (see also [START_REF] Emako | Well-balanced and asymptotic preserving schemes for kinetic models[END_REF])

ε∂ t f + v∂ x f = 1 ε V T ε (t, x, v ′ )f (t, x, v ′ ) dv ′ 2 -T ε (t, x, v)f (t, x, v) . (5.56)
The tumbling rate T ε describes the response to variations of chemical concentration along the path of bacteria. Among several choices, [START_REF] Calvez | Traveling chemotactic aggregates at mesoscopic scale and bi-stability[END_REF][START_REF] Othmer | The diffusion limit of transport equations II: Chemotaxis equations[END_REF], we choose T ε (t, x, v) = 1 + εφ v∂ x S(t, x) , φ an odd function.

(5.57)

In applications, the quantity S models the concentration of the chemoattractant which is released by bacteria themselves. It is then computed thanks to an elliptic/parabolic equation depending on the density of bacteria. Since we only focus on the diffusive limit of the kinetic system, we will consider that S is given, which boils down, from a numerical point of view, to treat explicitely in time the equation for S. Formally, the limit ε → 0 may be obtained easily by performing a Hilbert expansion, f = f 0 + εf 1 + . . ., equating each term in power of ε in (5.56),

f 0 = 1 2 ρ 0 , ρ 0 = V f 0 (v) dv; f 1 - 1 2 V f 1 (v ′ ) dv ′ = 1 2 V φ(v ′ ∂ x S)f 0 (v ′ ) dv ′ -φ(v∂ x S)f 0 -v∂ x f 0 .
Since φ is odd, V = (-1, 1) is symmetric and f 0 is independent of v, the first term of the right hand side vanishes. By conservation, we have

∂ t V f 0 (v) dv + ∂ x V vf 1 (v) dv = 0.
Along with the expressions of f 0 and f 1 ,

∂ t ρ 0 -∂ x 1 3 ∂ x ρ 0 + Eρ 0 = 0, E = 1 2 V vφ(v∂ x S) dv.
(5.58)

Sharfetter-Gummel scheme

The Il'in/Sharfetter-Gummel scheme (1.4)-(1.5) for this latter equation reads

ρ n+1 j = ρ n j + ∆t ∆x J n j-1 2 -J n j+ 1 2 , J n j-1/2 = E j-1 2 e 3E j-1 2 ∆x ρ n j -ρ n j-1 1 -e 3E j-1 2 ∆x , (5.59)
where E j-1 2 is a discretization of E at each interface of the mesh.

Remark 5.1 Clearly, if φ ≡ 0, the former model (4.29) is recovered out of (5.56), along with its limit (4.30), being a particular case of (5.58). Accordingly, (4.55) appears as a restriction of (5.59) when E j-1 2 ≡ 0. However, the situation φ ≡ 0 gives rise to sufficiently strong peculiarities so that we choose, in this paper, to neatly distinguish between both cases.

As in the previous Section, V = (-1, 1), the set {ω k , v k } k is assumed to verify (4.31).

Asymptotic expansion of eigenvalues

Let us focus on the eigenvalues of the discrete problem, computed thanks to the condition in (4.33), which, for ν =

1 2 K k=-K ω k δ(v -v k ), reads 1 = 1 2 K k=-K ω k T ε (v k ) T ε (v k ) -λ ε v k . ( 5.60)
Clearly, λ = 0 is a solution. Thanks to (4.31), we deduce from the latter equality that nonzero eigenvalues verify

N k=-N ω k Tε(v k ) v k -λ ε = 0, for λ ε = 0.
By studying the variations of the left-hand side with respect to λ ε , one deduces the existence of exactly 2K -1 distinct solutions which are interlaced in the following way

T ε (v -1 ) v -1 < λ ε -K+1 < T ε (v -2 ) v -2 < . . . < T ε (v -K ) v -K < λ ε 0 < T ε (v K ) v K < λ ε 1 < T ε (v K-1 ) v K-1 < λ ε 2 < . . . < λ ε K-1 < T ε (v 1 ) v 1 .
The sign of λ ε 0 is given by the sign of

K k=-K ω k v k T ε (v k ) .
In the following, we always assume that λ ε 0 < 0, the opposite case can be treated in the same way. Then, the vectors of negative/positive eigenvalues are,

λ ε -= (λ ε -K+1 , . . . , λ ε -1 ) ⊤ , λ ε + = (λ ε 1 , . . . , λ ε K-1 ) ⊤ . Lemma 5.2 When ε → 0, we have λ ℓ = λ 0 ℓ + ελ 1 ℓ + o(ε)
where

λ 0 0 = 0, λ 1 0 = 3 K k=1 ω k v k φ(v k ∂ x S).
For ℓ = 0, λ 0 ℓ are the symmetric (λ 0 -ℓ = -λ 0 ℓ ) eigenvalues of

1 = 1 2 K k=-K ω k 1 -λ 0 ℓ v k , and
λ 1 ℓ K k=-K ω k v k (1 -λ 0 ℓ v k ) 2 = λ 0 ℓ K k=-K ω k v k φ(v k ∂ x S) (1 -λ 0 ℓ v k ) 2 .
We denote by λ 0 = (λ 0 1 , . . . , λ 0 K-1 ) ⊤ the vector of positive eigenvalues at the limit ε → 0.

Proof. Letting ε → 0 in (5.60), eigenvalues λ 0 are solutions of

1 = 1 2 K k=-K ω k 1 -λ 0 v k .
Symmetry of the set {ω k , v k } k implies that if λ 0 is solution, then -λ 0 is, too. As a consequence λ 0 + = -λ 0 -, and λ 0 0 = 0. Assuming an asymptotic expansion λ ε k = λ 0 k + ελ 1 k + o(ε), and expanding the relation (5.60) with (5.57), we get

1 = 1 2 K k=-K ω k 1 -λ 0 v k -εv k (λ 1 -φ(v k ∂ x S)λ 0 ) + o(ε) = 1 2 K k=-K ω k 1 -λ 0 v k 1 + ε v k (λ 1 -φ(v k ∂ x S)λ 0 ) 1 -λ 0 v k + o(ε) .
Thus, for ℓ = 0,

λ 1 ℓ K k=-K ω k v k (1 -λ 0 ℓ v k ) 2 = λ 0 ℓ K k=-K ω k v k φ(v k ∂ x S) (1 -λ 0 ℓ v k ) 2 .
For ℓ = 0, this relation gives λ 0 0 = 0, forcing us to go at the second order in ε to compute λ 1 0 : postulating that

λ ε 0 = ελ 1 0 + ε 2 λ 2 0 + . . ., it comes 1 = 1 2 K k=-K ω k 1 -εv k (λ 1 0 + ελ 2 0 )/(1 + εφ(v k ∂ x S)) = 1 2 K k=-K ω k 1 + εv k λ 1 0 + ελ 2 0 1 + εφ(v k ∂ x S) + ε 2 v 2 k (λ 1 0 ) 2 + o(ε 2 ) . We get 0 = K k=-K ω k v k λ 1 0 + ε K k=-K ω k -v k φ(v k ∂ x S)λ 1 0 + v k λ 2 0 + v 2 k (λ 1 0 ) 2 + o(ε).
By symmetry of the set {ω k , v k }, we have K k=-K ω k v k = 0. Then, assumptions (4.31) yield two solutions, among which we discard the null one, it gives

λ 1 0 = 3 2 K k=-K ω k v k φ(v k ∂ x S). 20
Being φ an odd function and by symmetry of {ω k , v k } k , the claim is proved.

Corresponding scattering S-matrix

A general stationary solution reads, (see the separation of variables in [START_REF] Calvez | Confinement by biased velocity jumps: aggregation of Escherichia coli[END_REF])

f ε (x, v) = K-1 ℓ=1 a ℓ e -λ ε ℓ x/ε T ε (v) -λ ε ℓ v + a T ε (v) + -1 ℓ=-K+1 a ℓ e -λ ε ℓ (x-∆x)/ε T ε (v) -λ ε ℓ v + a 0 e -λ ε 0 x ε T ε (v) -λ ε 0 v - 1 T ε (v)
.

The spectral component of the "zero-eigenfunction" 1/T ε (v) was split between a and a 0 . Like for radiative transfer (see §4.2), where φ ≡ 0, the S-matrix is,

S ε = N ε (N ε ) -1 , (dependent of j -1 2 ) (5.61)
where the j -1 2 index was dropped for the sake of simplicity of the scripture,

N ε =   1 Tε(V)-V⊗λ + 1 Tε(V) e λ -∆x/ε Tε(V)-V⊗λ - 1 Tε(V)-λ ε 0 V -1 Tε(V) e -λ + ∆x/ε Tε(-V)+V⊗λ + 1 Tε(-V) 1 Tε(-V)+V⊗λ - e -λ ε 0 ∆x/ε Tε(-V)+λ ε 0 V -1 Tε(-V)   , N ε =   e -λ + ∆x/ε Tε(V)-V⊗λ + 1 Tε(V) 1 Tε(V)-V⊗λ - e -λ ε 0 ∆x/ε Tε(V)-λ ε 0 V -1 Tε(V) 1 Tε(-V)+V⊗λ + 1 Tε(-V) e λ -∆x/ε Tε(-V)+V⊗λ - 1 Tε(-V)+λ ε 0 V -1 Tε(-V)   .
The main differences with radiative transfer, as in §4, are that 0 is a simple eigenvalue but becomes a double eigenvalue when ε → 0, and that all above quantities depend on the index j, which lead to more intricate computations.

Decomposition of the scattering matrix

We want to compute the expansion (2.11) for this scattering matrix.

Lemma 5.3

The scattering matrix for the kinetic model for chemotaxis, defined in (5.61) admits the following asymptotic expansion in ε,

S ε = 0 K I K -ζ 0 γ I K -ζ 0 γ 0 K + εB ε , B ε = 1 ε (A ε (N ε ) -1 -A 0 (N 0 ) -1 ), ( 5

.62)

where

A ε = N ε - 0 K I K I K 0 K N ε
, and where the matrices γ ∈ M K-1×K (R) and

ζ 0 ∈ M K×K-1 (R) satisfy, γ 1 1 -V ⊗ λ 0 = I K-1 , γ1 R K = 0 R K-1 , (5.63)
ζ 0 = 1 1 -V ⊗ λ 0 - 1 1 + V ⊗ λ 0 .
(5.64) Remark 5.4 Lemma 5.3 is the equivalent of Lemma 4.3 in the case of the radiative transfert equation. The existence of γ is also provided by Proposition A.1 (iii) in Appendix.

Proof. The index j -1/2 is again dropped since there is no possible confusion.

1. When ε → 0, being T ε given in (5.57) and using Lemma 5.2

1 T ε (V) -λ ε 0 V - 1 T ε (V) → 0, e -λ ε 0 ∆x/ε T ε (V) -λ ε 0 V - 1 T ε (V) → e -λ 1 0 ∆x -1.
so that S ε → S 0 , where S 0 = N 0 (N 0 ) -1 . Since λ 0 := λ 0 + = -λ 0 -,

N 0 = 1 1-V⊗λ 0 1 R K 0 K 0 K×K-1 1 R K 1 1-λ 0 ⊗V (e -λ 1 0 ∆x -1)1 R K , N 0 = 0 K×K-1 1 R K 1 1+V⊗λ 0 (e -λ 1 0 ∆x -1)1 R K 1 1+λ 0 ⊗V 1 R K 0 K .
2. Using γ defined in (5.63), we define also β ∈ R K such that

β ⊤ 1 1 -V ⊗ λ 0 = 0 ⊤ R K , β ⊤ 1 R K = 1.
(5.65)

Then, we may write

(N 0 ) -1 =      γ β ⊤ 0 K - γ 1 e -λ 1 0 ∆x -1 β ⊤ 0 K×K-1 1 R K γ β ⊤ γ 1 e -λ 1 0 ∆x -1 β ⊤      =      γ β ⊤ 0 K - 0 K-1×K 1 e -λ 1 0 ∆x -1 β ⊤ γ 1 e -λ 1 0 ∆x -1 β ⊤      .
By definition,

A ε = N ε - 0 K I K I K 0 K N ε , and so S ε = 0 K I K I K 0 K + A ε (N ε ) -1 .

Let us denote

ζ ε ± = ±1 T ε (V) -V ⊗ λ ε ± - ±1 T ε (-V) + V ⊗ λ ε ± ∈ M K×K-1 (R), ζ ε 0 = 1 T ε (V) -λ ε 0 V - 1 T ε (-V) + λ ε 0 V ∈ R K .
As ε → 0, we obtain the limit

ζ ε ± → ζ 0 = 1 1 -V ⊗ λ 0 - 1 1 + V ⊗ λ 0 ,
and

ζ ε 0 → 0.
From the expression,

A ε = e -λ + ∆x/ε ζ ε + 1 Tε(V) -1 Tε(-V) -ζ ε - e -λ ε 0 ∆x/ε ζ ε 0 -1 Tε(V) + 1 Tε(-V) -ζ ε + 1 Tε(-V) -1 Tε(V) e λ -∆x/ε ζ ε - -ζ ε 0 + 1 Tε(V) -1 Tε(-V)
.

(5.66) in the limit ε → 0, we get

A 0 =   0 K -ζ 0 0 R K -ζ 0 0 R K 0 K   , A 0 (N 0 ) -1 = 0 K -ζ 0 γ -ζ 0 γ 0 K . (5.67)
Thus we reach decomposition (5.62).

Emergence of an asymptotic scheme

We deduce the final scheme from (2.12), which mostly reads as (4.49),

f n+1 j (V) f n+1 j-1 (-V) + ∆t ε∆x V f n+1 j (V) -(I K -ζ 0 j-1 2 γ j-1 2 )f n+1 j (-V) f n+1 j-1 (-V) -(I K -ζ 0 j-1 2 γ j-1 2 )f n+1 j-1 (V) = f n j (V) f n j-1 (-V) + ∆t ∆x VB ε j-1 2 f n j-1 (V) f n j (-V)
.

(5.68)

Denoting B ε = B 1ε B 2ε B 3ε B 4ε
, we may rewrite (5.68) as

1 ε R ε j f n+1 j (V) f n+1 j (-V) = f n j (V) f n j (-V)) + ∆t ∆x V B 1ε j-1 2 f n j-1 (V) + B 2ε j-1 2 f n j (-V) B 3ε j+ 1 2 f n j (V) + B 4ε j+ 1 2 f n j+1 (-V) , (5.69)
where

R ε j = εI 2K + ∆t ∆x V I K ζ 0 j-1 2 γ j-1 2 -I K ζ 0 j+ 1 2 γ j+ 1 2 -I K I K . (5.70)
After inversion of the matrix R ε , by construction, this scheme satisfies both the well-balanced property and is asymptotic preserving.

Our main result for the Othmer-Alt model for chemotaxis (5.56) reads:

Theorem 5.5 The scheme (5.69) for the chemotaxis model (5.56) is well-balanced and uniformly accurate (AP) with respect to ε. Moreover, when ε → 0, the macroscopic density ρ n j :=

K k=-K ω k f n j (v k )
satisfies the Sharfetter-Gummel discretization (5.59), where E j-

1 2 := K k=1 ω k v k φ(v k ∂ x S j-1 2
).

The proof of this Theorem is done in the next subsection.

Consistency with the limit ε → 0

We first state the following technical Lemma (dropping the subscript j -1 2 since there is no possible confusion): Lemma 5.6 When ε goes to 0,

B ε = B 1ε B 2ε B 3ε B 4ε → B 0 = B 10 B 20 B 30 B 40 , where                            B 10 = -e -λ 1 0 ∆x e -λ 1 0 ∆x -1 (2 -ζ 0 γ)φ(V∂ x S)β ⊤ ; B 20 = (ζ 0 γ -1) φ(V∂xS)-V⊗λ 1 - (1-V⊗λ 0 ) 2 + φ(V∂xS)-V⊗λ 1 - (1+V⊗λ 0 ) 2 γ + 1 e -λ 1 0 ∆x -1 (2 -ζ 0 γ)φ(V∂ x S)β ⊤ ; B 30 = (ζ 0 γ -1) -φ(V∂xS)+V⊗λ 1 + (1-V⊗λ 0 ) 2 + -φ(V∂xS)+V⊗λ 1 + (1+V⊗λ 0 ) 2 γ + (2 -ζ 0 γ)φ(V∂ x S)β ⊤ + 1 e -λ 1 0 ∆x -1 (2 -ζ 0 γ)λ 1 0 Vβ ⊤ ; B 40 = - 1 e -λ 1 0 ∆x -1 (2 -ζ 0 γ)λ 1 0 Vβ ⊤ ;
with γ, ζ 0 , and β defined in (5.63), (5.64), and (5.65) respectively.

Proof.

1. First, rewrite 5.66) and (5.67), we have

B ε = 1 ε (A ε -A 0 )N -1 ε + 1 ε (A 0 -A 0 N -1 0 N ε )N -1 ε . From (
A ε -A 0 ε =   1 ε e -λ + ∆x/ε ζ ε + 1 ε 1 Tε(V) -1 Tε(-V) -1 ε (ζ ε --ζ 0 ) -1 ε (ζ ε + -ζ 0 ) 1 ε 1 Tε(-V) -1 Tε(V) 1 ε e λ -∆x/ε ζ ε - 1 ε e -λ ε 0 ∆x/ε ζ ε 0 -1 Tε(V) + 1 Tε(-V) 1 ε -ζ ε 0 + 1 Tε(V) -1 Tε(-V)   ,
where, from (5.57),

lim ε→0 1 ε 1 T ε (V) - 1 T ε (-V) = -2φ(V∂ x S). 2. Denoting λ 1 ± = lim ε→0 1 ε (λ ε ± -λ 0 ) (as in Lemma 5.2), we get lim ε→0 1 ε (ζ ε ± -ζ 0 ) = V ⊗ λ 1 ± -φ(V∂ x S) (1 ∓ V ⊗ λ 0 ) 2 - V ⊗ λ 1 ± -φ(V∂ x S) (1 ± V ⊗ λ 0 ) 2 . Also, 1 ε ζ ε 0 → -2φ(V∂ x S) + 2λ 1 0 V. Thus, when ε → 0, lim ε→0 A ε -A 0 ε = 0 K×K-1 -2φ(V∂ x S) φ(V∂xS)-V⊗λ 1 + (1-V⊗λ 0 ) 2 - φ(V∂xS)-V⊗λ 1 + (1+V⊗λ 0 ) 2 2φ(V∂ x S) φ(V∂xS)-V⊗λ 1 - (1+V⊗λ 0 ) 2 - φ(V∂xS)-V⊗λ 1 - (1-V⊗λ 0 ) 2 2φ(V∂ x S) 0 K×K-1 -2λ 1 0 V . 3. Thanks to the fact that γ 1 1-V⊗λ 0 = I K-1 , it comes 1 ε (A 0 -A 0 N -1 0 N ε ) =   ζ 0 γ ε e -λ + ∆x/ε Tε(-V)+V⊗λ + 1 Tε(-V) ζ 0 γ ε 1 Tε(V)-V⊗λ + - 1 1-V⊗λ 0 1 Tε(V) ζ 0 γ ε 1 Tε(-V)+V⊗λ -- 1 1-V⊗λ 0 e -λ ε 0 ∆x/ε Tε(-V)+λ ε 0 V -1 Tε(-V) ζ 0 γ ε e λ -∆x/ε Tε(V)-V⊗λ - 1 Tε(V)-λ ε 0 V -1 Tε(V)   . Since γ1 R K = 0 R K-1 , we have γ 1 εT ε (V) = γ ε 1 T ε (V) -1 R K = γ -φ(V∂ x S) 1 + εφ(V∂ x S) -→ ε→0 -γφ(V∂ x S).
We have also

γ ε 1 T ε (V) -V ⊗ λ + - 1 1 -V ⊗ λ 0 -→ ε→0 γ -φ(V∂ x S) + V ⊗ λ 1 + (1 -λ 0 ⊗ V) 2 .
Then,

lim ε→0 1 ε (A 0 -A 0 N -1 0 N ε ) =   ζ 0 γ (0 K×K-1 φ(V∂ x S)) ζ 0 γ φ(V∂xS)-V⊗λ 1 - (1-λ 0 ⊗V) 2 -φ(V∂ x S) ζ 0 γ -φ(V∂xS)+V⊗λ 1 + (1-λ 0 ⊗V) 2 -φ(V∂ x S) ζ 0 γ 0 K×K-1 λ 1 0 V   .
4. Gathering these computations, we reach

B 0 := lim ε→0 B ε = lim ε→0 A ε -A 0 ε (N 0 ) -1 + lim ε→0 1 ε (A 0 -A 0 (N 0 ) -1 N ε )(N 0 ) -1 = B 10 B 20 B 30 B 40 ,
where B 10 , B 20 , B 30 , and B 40 are expressed in Lemma 5.6.

Proof of Theorem 5.5 We proceed as in Section 4.4. First, we have obviously from (5.70)

R ε j = εI 2K + R 0 j where R 0 j := ∆t ∆x V I K ζ 0 j-1 2 γ j-1 2 -I K ζ 0 j+ 1 2 γ j+ 1 2 -I K I K
Assuming that f admits an Hilbert expansion f = f 0 + εf 1 + . . ., we get, by injecting this expansion into (5.69) and identifying the term in power of ε,

R 0 j {f 0 } n+1 j (V) {f 0 } n+1 j (-V) = 0, (5.71) and R 0 j {f 1 } n+1 j (V) {f 1 } n+1 j (-V) = 1 2 {f 0 } n j -{f 0 } n+1 j {f 0 } n j -{f 0 } n+1 j (5.72) + ∆t 2∆x V B 10 j-1 2 {f 0 } n j-1 + B 20 j-1 2 {f 0 } n j B 30
j+ 1 2 {f 0 } n j + B 40 j+ 1 2 {f 0 } n j+1 .
Since from (4.36) we have that ∀ℓ, K k=1 ω k v k (ζ 0 γ) kℓ = 0, we may apply Lemma A.2 and deduce

• Ker(R 0 j ) = Span{1 R K }, • Im(R 0 j ) = {Z = (Z 1 Z 2 ) ⊤ , Z i ∈ R K such that K k=1 ω k (Z 1k + Z 2k ) = 0}. Thus equation (5.71) implies {f 0 } n+1 j (±V) = ρ n+1 j 2 1 R K .
Injecting into (5.72), we get

R 0 j {f 1 } n+1 j (V) {f 1 } n+1 j (-V) = 1 2 (ρ n j -ρ n+1 j )1 R K (ρ n j -ρ n+1 j )1 R K (5.73) + ∆t 2∆x V B 10 j-1 2 1 R K ρ n j-1 + B 20 j-1 2 1 R K ρ n j B 30
j+ 1 2 1 R K ρ n j + B 40 j+ 1 2 1 R K ρ n j+1 .
Thanks to the relations,

γ1 R K = 0 R K-1 . β ⊤ 1 R K = 1,
we deduce from the expression of B 0 in Lemma 5.6,

B 10 j-1 2 1 R K = - e -λ 1 0,j-1 2 ∆x e -λ 1 0,j-1 2 ∆x -1 (2 -ζ 0 γ)φ(V∂ x S j-1 2
);

B 20 j-1 2 1 R K = 1 e -λ 1 0,j-1 2 ∆x -1 (2 -ζ 0 γ)φ(V∂ x S j-1 2
);

B 30 j+ 1 2 1 R K = (2 -ζ 0 γ)φ(V∂ x S j+ 1 2 ) + λ 1 0,j+ 1 2 e -λ 1 0,j+ 1 2 ∆x -1 (2 -ζ 0 γ)V; B 40
j+ 1 2 1 R K = -λ 1 0,j+ 1 2 e -λ 1 0,j+ 1 2 ∆x -1 (2 -ζ 0 γ)V.
The solvability condition for (5.73) is that its right hand side belongs to Im(R 0 j ). Then, we multiply each line by ω k and add over k to obtain

0 =ρ n j -ρ n+1 j + ∆t ∆x   -e -λ 1 0,j-1 2 ∆x ρ n j-1 + ρ n j e -λ 1 0,j-1 2 ∆x -1 K k=1 ω k v k φ(v k ∂ x S j-1 2 ) +ρ n j K k=1 ω k v k φ(v k ∂ x S j+ 1 2 ) + λ 1 0,j+ 1 2 (ρ n j -ρ n j+1 ) e -λ 1 0,j+ 1 2 ∆x -1 K k=1 ω k v 2 k   ,
where we have used also ∀ℓ, K k=1 ω k v k (ζ 0 γ) kℓ = 0. From Lemma 5.2 and assumptions (4.31), it comes 1 3

λ 1 0 = K k=1 ω k v k φ(v k ∂ x S)
.

Denoting E j-1 2 := K k=1 ω k v k φ(v k ∂ x S j-1 2
), the above scheme rewrites as

ρ n+1 j = ρ n j + ∆t ∆x   -E j-1 2 e -3E j-1 2 ∆x ρ n j-1 -ρ n j e -3E j-1 2 ∆x -1 + E j+ 1 2 e -3E j+ 1 2 ∆x ρ n j -ρ n j+1 e -3E j+ 1 2 ∆x -1   ,
in which we recognize the Sharfetter-Gummel scheme (5.59).

6 Diffusive limit of Vlasov-Fokker-Planck kinetic equations

Presentation of the continuous model

This kinetic model reads, see [START_REF] Nieto | High field limit for the Vlasov-Poisson-Fokker-Planck system[END_REF][START_REF] Poupaud | Parabolic limit and stability of the Vlasov-Poisson-Fokker-Planck system[END_REF][START_REF] Wollman | Numerical approximation of the Vlasov-Poisson-Fokker-Planck system in one dimension[END_REF], in parabolic scaling,

ε∂ t f ε + v∂ x f ε + E • ∂ v f ε = 1 ε ∂ v (vf ε + κ∂ v f ε ) , 0 < ε ≪ 1. (6.74)
The electric field is denoted E and may derive from a potential φ through the relation

E = ±∂ x φ.
When such a potential φ depends self-consistently on the macroscopic density of electrons (through the Poisson equation), one speaks about the Vlasov-Poisson-Fokker-Planck system. To keep the exposition simple, we shall assume that a steady electric field E(x) is given. When ε → 0, the kinetic distribution in (6.74) relaxes to a Maxwellian,

f ε → ρ 0 1 √ 2πκ e -v 2 /2κ ,
where the macroscopic density ρ 0 solves a drift-diffusion (continuity) equation,

∂ t ρ 0 + κ∂ x E κ ρ 0 -∂ x ρ 0 = 0. ( 6.75)
The Sharfetter-Gummel scheme associated to this system reads

ρ n+1 j = ρ n j + ∆t ∆x J n j-1 2 -J n j+ 1 2 , (6.76)
J n j-1/2 = E j-1 2 ρ n j-1 -e -E j-1 2 ∆x/κ ρ n j 1 -e -E j-1 2 ∆x/κ .

Spectral decomposition of stationary solutions

Consider the Fokker-Planck stationary problem with inflow boundaries,

v∂ x f = 1 ε ∂ v (v -εE)f + κ∂ v f , ε, κ > 0. (6.77)
A convenient "separated variables" ansatz reads now, (see e.g. [START_REF] Burschka | The Kinetic Boundary Layer for the Fokker-Planck Equation with Absorbing Boundary[END_REF])

f (x, v) = exp(-λx -µv)ψ λ (v),
so that one recovers a standard Sturm-Liouville eigenvalue problem (see [START_REF] Beals | Partial-range completeness and existence of solutions to two-way diffusion equations[END_REF] and [START_REF] Gosse | Computing Qualitatively Correct Approximations of Balance Laws[END_REF]Chapter 12.3]) with a discrete spectrum. The null eigenvalue λ = 0 is double, its two associated non-damped modes are denoted Ψ ε ±0 ("diffusion solutions" in [START_REF] Fisch | Separating variables in two-way diffusion equations[END_REF]) among which appears a space-homogeneous mode:

Ψ ε ±0 (x, v) = exp(- µ ε ±0 x ε )ψ ε ±ℓ (v),
where

• when E > 0,

µ ε 0 = 0; ψ ε 0 (v) = exp(- (v -|εE|) 2 2κ
);

µ ε -0 = - εE κ ; ψ ε -0 (v) = exp(- |εE| 2 2κ ) exp(- v 2 2κ
);

• when E < 0,

µ ε 0 = - εE κ ; ψ ε 0 (v) = exp(- |εE| 2 2κ ) exp(- v 2 2κ
);

µ ε -0 = 0; ψ ε -0 (v) = exp(- (v + |εE|) 2 2κ
).

Other eigenfunctions Ψ ε ±ℓ , for ℓ ∈ N * , are explicitly given in [24, p. 251]:

Ψ ε ±ℓ (x, v) = exp(- µ ε ±ℓ x ε ) ψ ε ±ℓ (v) , (6.78)
ψ ε ±ℓ (v) = exp(-µ ε ±ℓ v)H ℓ (ṽ ε ±ℓ ) exp(-(ṽ ε ±ℓ ) 2 ) , (6.79)
where H ℓ is the ℓth Hermite polynomial and

µ ε ±ℓ = -εE ± (εE) 2 + 4κℓ 2κ , ṽε ±ℓ = v -2µ ε ±ℓ κ -εE √ 2κ .
Hermite's polynomials are such that deg(H ℓ ) = ℓ and H ℓ (-X) = (-1) ℓ H ℓ (X). They are orthogonal with respect to a strongly-growing (indefinite) weight:

R v ψ ε ±k (v)ψ ε ±ℓ (v) exp (v -εE) 2 2κ dv = 0, if k = ℓ. (6.80)
A spectral decomposition follows, for smooth enough functions f (x, v) [START_REF] Beals | Half-range completeness for the Fokker-Planck equation[END_REF][START_REF] Cercignani | Half-Range Completeness for the Fokker-Planck Equation with an External Force[END_REF],

f (x, v) = α +0 Ψ ε 0 (x, v) + α -0 Ψ ε -0 (x, v) + ℓ≥1 α ℓ Ψ ε ℓ (x, v) + α -ℓ Ψ ε -ℓ (x, v) . (6.81)
By inserting the space-homogeneous eigenfunction associated to k = 0, one sees from (6.80) that none in all the set of Ψ ε ±ℓ , ℓ > 0, can carry any net macroscopic flux:

∀ℓ = 0, R vψ ε ±ℓ (v)ψ ε ±0 (v) exp (v -εE) 2 2κ dv = R vψ ε ±ℓ (v)dv = 0. (6.82)
For future use, we compute the limit as ε → 0 of the above expressions. We obtain straightforwardly, as ε → 0,

ψ 0 ±0 (v) = exp - v 2 2κ ; (6.83a) ψ 0 ±ℓ (v) = H ℓ v ∓ 2 √ ℓκ √ 2κ exp - v 2 2κ ± v ℓ κ -2ℓ , for ℓ ∈ N * . (6.83b)
We deduce from the symmetry of Hermite polynomials that we have the identity ψ 0 -ℓ (-v) = (-1) ℓ ψ 0 ℓ (v). These expressions do not depend on E ∈ R.

Assumptions on the set of discrete velocities and weights.

As mentioned in the introduction, and contrary to the Case's functions, being exponential polynomials, solutions of the stationary Vlasov-Fokker-Planck equation, do not constitute a Chebyshev T -system on (0, +∞). Indeed, if we assume that the family {ψ 0 0 , . . . , ψ 0 K-1 } is endowed with the Haar property, then any non-zero linear combination of these functions will have no more than K -1 roots on (0, +∞); otherwise, denoting v 0 , . . . , v K-1 these roots, det(ψ 0 ℓ (v k )) k,ℓ = 0 (see Proposition 2.3). However, it is not difficult to find ad-hoc coefficients a 0 , . . . , a K-1 , for which the function

R + ∋ v → a 0 ψ 0 0 (v) + K-1 ℓ=1 a ℓ ψ 0 ℓ (v)
admits more than K -1 roots on (0, +∞). For instance,

• for K = 2, the function (after simplification by exp(-

v 2 /2κ)), v → 3 2 + v -2 √ 2 e v
has two positive roots, approximately given by 0.1216 and 1.5495;

• For K = 3, the function v → -2.75 + 0.2e v (v -2) + e √ 2v v 2 2 -2 √ 2v + 3
admits three positive roots, approximately given by 0.132, 0.796 and 4.192.

The number of real roots for exponential polynomials admits the Pólya-Szegö estimate as an upper bound, see Appendix 3.1. However, exponential monomials do satisfy the Haar property on (0, +∞) as shown in Theorem A.1. Accordingly, we must prescribe some assumptions on the set of discrete velocities. More precisely, we assume that v 1 , . . . , v K are chosen such that,

given K ∈ N, det ψ 0 0 (V) ψ 0 1 (V) . . . ψ 0 K-1 (V) = 0, (6.84)
and such that the family

ψ 0 0 (V) ψ 0 0 (-V) , ψ 0 1 (V) ψ 0 1 (-V) , • • • , ψ 0 K-1 (V) ψ 0 K-1 (-V) , ψ 0 1 (-V) ψ 0 1 (V) , • • • , ψ 0 K-1 (-V) ψ 0 K-1 (V) is linearly independent in R 2K . (6.85)
We assume moreover that the corresponding weights ω k are such that the orthogonality relation (6.82) holds true at the discrete level, i.e.

∀ℓ = 1, . . . , K -1, K k=1 ω k v k (ψ 0 ±ℓ (v k ) -ψ 0 ±ℓ (-v k )) = 0. (6.86) Thanks to the relation ψ 0 -ℓ (-v k ) = (-1) ℓ ψ 0 ℓ (v k )
, it suffices to get the above identity only for positive ℓ.

Finally we mention the following Lemma which shows that actually assumptions (6.84) and (6.85) hold at the continuous level. Lemma 6.1 For any N ∈ N, the family

(ψ 0 0 , ψ 0 1 , • • • , ψ 0 N -1 )
given by (6.79) is linearly independent on v ∈ (0, +∞), and also on v ∈ R.

Proof. Given N , we have to show the linear independence on R and on (0, +∞). The result on R is a direct consequence of the orthogonality relation (6.80). Thus we are left to prove the result on (0, +∞). After expanding the Gaussian term in (6.79) and normalizing coefficients, any linear combination with ε = 0 rewrites, ∀v ∈ (0, +∞),

N -1 i=0 λ i exp(µ i v)H i (ṽ i ) = 0. Since H N -1 = 0 for v big enough, it suffices to let v → +∞ in λ N -1 + N -2 i=0 λ i exp((µ i -µ N -1 ) v) H i (ṽ i ) H N -1 (ṽ N -1 )
= 0, in order to get λ N -1 = 0 because µ iµ N -1 ≤ C < 0. Successive coefficients vanish for the same reason. Hence for any N , the family (ψ 0 i ) i<N is linearly independent on R + * .

Corresponding scattering S-matrix

The stationary problem with incoming boundary condition reads

v∂ x f = 1 ε ∂ v (v -εE j- 1 2
)f + κ∂ v f , on (0, ∆x), (6.87)

f (0, v) = f n j-1 (v), f (∆x, -v) = f n j (-v). (6.88)
Based on (6.81) and (6.78), we seek spectral approximations obtained by truncating (6.81) to the first 2K modes,

f (x, v) = α +0 ψ ε 0 (v)e -µ ε 0 x/ε + α -0 ψ ε -0 (v)e -µ ε -0 x/ε + K-1 ℓ=1 α ±ℓ ψ ε ±ℓ (v)e -µ ε ±ℓ x/ε . (6.89)
Coefficients α ±ℓ in this full-range expansion follow from incoming boundary conditions (6.88), according to a linear system: for k = 1, . . . , K,

f (0, v k ) = f n j-1 (v k ) =α 0 ψ ε 0 (v k ) + α -0 ψ ε -0 (v j ) + K-1 ℓ=1 α ±ℓ ψ ε ±ℓ (v k ), f (∆x, -v k ) = f n i+1 (-v k ) =α 0 ψ ε 0 (-v k )e -µ 0 ε ∆x + α -0 ψ ε -0 (-v k )e - µ ε -0 ε ∆x + K-1 ℓ=1 α ±ℓ ψ ε ±ℓ (-v k )e - µ ε ±ℓ ε ∆x .
The K × (K -1) matrix of eigenvectors associated to nonzero eigenvalues is

ψ ε ± (V) = ψ ε ±1 (V) . . . ψ ε ±K-1 (V) , so that, α + α - = (M ε ) -1 f n j-1 (V) f n j (-V)
, where M ε is the 2K × 2K matrix defined by

M ε = ψ ε + (V) ψ ε 0 (V) ψ ε -(V) ψ ε -0 (V) ψ ε + (-V)e -µ ε + ∆x/ε ψ ε 0 (-V)e -µ ε 0 ∆x/ε ψ ε -(-V)e -µ ε -∆x/ε ψ ε -0 (-V)e -µ ε -0 ∆x/ε .
Outgoing values of f follow thanks to (6.89),

f (∆x, V) = α 0 ψ ε 0 (V)e -µ ε 0 ∆x/ε + α -0 ψ ε -0 (V)e -µ ε -0 ∆x/ε + K-1 ℓ=1 α ±ℓ ψ ε ±ℓ (V)e -µ ε ±ℓ ∆x/ε , f (0, -V) = α 0 ψ ε 0 (-V) + α -0 ψ ε -0 (-V) + K-1 ℓ=1 α ±ℓ ψ ε ±ℓ (-V).
Then, denoting

M ε = ψ ε + (V)e -µ ε + ∆x/ε ψ ε 0 (V)e -µ ε 0 ∆x/ε ψ ε -(V)e -µ ε -∆x/ε ψ ε -0 (V)e -µ ε -0 ∆x/ε ψ ε + (-V) ψ ε 0 (-V) ψ ε -(-V) ψ ε -0 (-V) , we get f (∆x, V) f (0, -V) = M ε (M ε ) -1 f n j-1 (V) f n j (-V)
.

Accordingly, the scattering matrix is defined by,

S ε j-1/2 = M ε (M ε ) -1 .
(dependent on j -1 2 )

Decomposition of the scattering matrix

We perform the computation in the case E > 0. Computations in the case E < 0 are similar. We start from the expression,

S ε = M ε (M ε ) -1 = N ε (N ε ) -1 ,
where (expressions of µ ε ±0 are used),

N ε = ψ ε + (V) ψ ε 0 (V) ψ ε -(V)e µ ε -∆x/ε ψ ε -0 (V)e -E∆x/κ ψ ε + (-V)e -µ ε + ∆x/ε ψ ε 0 (-V) ψ ε -(-V) ψ ε -0 (-V) , N ε = ψ ε + (V)e -µ ε + ∆x/ε ψ ε 0 (V) ψ ε -(V) ψ ε -0 (V) ψ ε + (-V) ψ ε 0 (-V) ψ ε -(-V)e µ ε -∆x/ε ψ ε -0 (-V)e -E∆x/κ .
Moreover, the product N ε (N ε ) -1 is invariant if, in both matrices N ε and N ε , we subtract to the last column e -E∆x/κ times the Kth column, i.e.

S ε = N ε (N ε ) -1 , (6.90)
with

N ε = ψ ε + (V) ψ ε 0 (V) ψ ε -(V)e µ ε -∆x/ε (ψ ε -0 (V) -ψ ε 0 (V))e -E∆x/κ ψ ε + (-V)e -µ ε + ∆x/ε ψ ε 0 (-V) ψ ε -(-V) ψ ε -0 (-V) -ψ ε 0 (-V)e -E∆x/κ , N ε = ψ ε + (V)e -µ ε + ∆x/ε ψ ε 0 (V) ψ ε -(V) ψ ε -0 (V) -ψ ε 0 (V)e -E∆x/κ ψ ε + (-V) ψ ε 0 (-V) ψ ε -(-V)e µ ε -∆x/ε (ψ ε -0 (-V) -ψ ε 0 (-V))e -E∆x/κ .
Noticing that when ε → 0,

ψ ε ±0 (V) → exp(-V 2 2κ
), we may pass to the limit ε → 0 in the latter matrices and get S ε → S 0 = N 0 ( N 0 ) -1 with

N 0 = ψ 0 + (V) exp(-V 2 2κ ) 0 K 0 K×K-1 exp(-V 2 2κ ) ψ 0 -(-V) (1 -e -E∆x/κ ) exp(-V 2 2κ ) , (6.91)
N 0 = 0 K×K-1 exp(-V 2 2κ ) ψ 0 -(V) (1 -e -E∆x/κ ) exp(-V 2 2κ ) ψ 0 + (-V) exp(-V 2 2κ ) 0 K . ( 6
.92) Lemma 6.2 Under assumption (6.84) on the set of discrete velocities, the S-matrix (6.90) admits the asymptotic expansion in ε,

S ε = 0 K I K -ζ 0 + γ + I K -ζ 0 + γ + 0 K + εB ε , B ε = 1 ε ( N ε (N ε ) -1 -N 0 (N 0 ) -1 ), (6.93)
where the matrices γ + ∈ M K-1×K (R) and

ζ 0 + ∈ M K×K-1 (R) satisfy, γ + =    γ ⊤ 1 . . . γ ⊤ K-1    where γ ℓ ∈ R K , γ ⊤ ℓ ψ 0 k (V) = δ kℓ , γ ⊤ ℓ exp(- V 2 2κ ) = 0; (6.94) ζ ε ± = ζ ε ±1 . . . ζ ε ±(K-1) , ζ ε ±ℓ = ψ ε ±ℓ (V) -ψ ε ±ℓ (-V) ∈ R K . ( 6
.95) Remark 6.3 Existence of γ + is guaranteed by assumption (6.84). Notice that since the limits ψ 0 ±ℓ in (6.83) do not depend on E, the first term in the decomposition of the scattering matrix is independent on E, hence on j.

Proof.

1. Let vector β ∈ R K be such that

β ⊤ ψ 0 ℓ (V) = 0; β ⊤ exp(- V 2 2κ ) = 1. (6.96)
Moreover, since

γ + =    γ ⊤ 1 . . . γ ⊤ K-1    , it comes that γ + β ⊤ = ψ 0 + (V) exp(-V 2 2κ ) -1 , (6.97)
where the latter matrix is invertible thanks to assumption (6.84). It brings also that,

K-1 ℓ=1 ψ 0 ℓ (V)γ ⊤ ℓ + exp(- V 2 2κ )β ⊤ = I K . (6.98) 2.
Similarly, as ψ 0 -ℓ (V) = (-1) ℓ ψ 0 ℓ (-V), then denoting

γ -=      -γ ⊤ 1 γ ⊤ 2 . . . (-1) K-1 γ ⊤ K-1      , it comes that, γ - 1 1-e -E∆x/κ β ⊤ = ψ 0 -(-V) (1 -e -E∆x/κ ) exp(-V 2 2κ ) -1
.

(6.99)

3. We deduce from both (6.97) and (6.99) that,

(N 0 ) -1 =     γ + β ⊤ 0 K 0 K-1×K - 1 1-e -E∆x/κ β ⊤ γ - 1 1-e -E∆x/κ β ⊤     ,
and so, with (6.91)-(6.92),

S ε -→ ε→0 S 0 := N 0 (N 0 ) -1 ,
has the form (6.93), with

N 0 (N 0 ) -1 =       0 K K-1 ℓ=1 ψ 0 ℓ (-V)γ ⊤ ℓ + exp(- V 2 2κ )β ⊤ K-1 ℓ=1 ψ 0 ℓ (-V)γ ⊤ ℓ + exp(- V 2 2κ )β ⊤ 0 K       .
Using (6.98) allows to complete the proof. where

B 10 = 1 1 -e -E∆x κ 2EV κ exp(- V 2 2κ )β ⊤ - 1 1 -e -E∆x κ ζ 0 + γ + EV κ - E 2 κ exp(- V 2 2κ )β ⊤ ; B 20 = dζ ε - dε |ε=0 + ζ 0 + γ + dψ ε -(-V) dε |ε=0 γ -- e -E∆x κ 1 -e -E∆x κ 2EV κ exp(- V 2 2κ )β ⊤ - 1 1 -e -E∆x κ ζ 0 + γ + E 2 κ + 2EV κ e -E∆x κ exp(- V 2 2κ )β ⊤ ; B 30 = ζ 0 + γ + dψ ε -(V) dε |ε=0 - dζ ε + dε |ε=0 γ + - 1 1 -e -E∆x κ 2EV κ exp(- V 2 2κ )β ⊤ + 1 1 -e -E∆x κ ζ 0 + γ + EV κ + E 2 κ exp(- V 2 2κ )β ⊤ ; B 40 = e -E∆x κ 1 -e -E∆x κ 2EV κ -ζ 0 + γ + E 2 κ + EV κ exp(- V 2 2κ )β ⊤ . Proof. Rewriting B ε = 1 ε ( N ε -N 0 (N 0 ) -1 N ε )(N ε ) -1 , it comes, B ε N ε = 1 ε ( N ε -N 0 (N 0 ) -1 N ε ) = 1 ε A 1ε A 2ε A 3ε A 4ε , with A 1ε = (ζ ε + ζ 0 + γ + ψ ε + (-V))e - µ ε + ∆x ε ζ ε 0 + ζ 0 + γ + ψ ε 0 (-V) ; A 2ε = ζ ε -+ ζ 0 + γ + ψ ε -(-V) ζ ε -0 + ζ 0 + γ + ψ ε -0 (-V) -(ζ ε 0 + ζ 0 + γ + ψ ε 0 (-V))e -E∆x κ ; A 3ε = ζ 0 + γ + ψ ε + (V) -ζ ε + ζ 0 + γ + ψ ε 0 (V) -ζ ε 0 ; A 4ε = (ζ 0 + γ + ψ ε -(V) -ζ ε -)e µ ε -∆x ε (ζ ε 0 -ζ ε -0 + ζ 0 + γ + (ψ ε -0 (V) -ψ ε 0 (V)))e -E∆x κ
.

When E > 0, ζ ε -0 = 0, whereas if E < 0, ζ ε 0 = 0. Accordingly, for E > 0, lim ε→0 1 ε ψ ε 0 (V) -exp(- V 2 2κ ) = dψ ε 0 (V) dε |ε=0 = EV κ exp(- V 2 2κ ); lim ε→0 1 ε ψ ε -0 (V) -exp(- V 2 2κ ) = - E 2 κ exp(- V 2 2κ ); lim ε→0 1 ε ζ ε 0 = 2EV κ exp(- V 2 2κ ),
and the expression (6.100) follows by doing B 0 = lim ε→0

1 ε A ε (N 0 ) -1 .

Uniform accuracy with respect to ε of the final scheme

We deduce the final scheme from (2.12), which reads like (5.68),

f n+1 j (V) f n+1 j-1 (-V) + ∆t ε∆x V f n+1 j (V) -(I K -ζ 0 + γ + )f n+1 j (-V) f n+1 j-1 (-V) -(I K -ζ 0 + γ + )f n+1 j-1 (V) = f n j (V) f n j-1 (-V) + ∆t ∆x VB ε j-1 2 f n j-1 (V) f n j (-V) . (6.101) Denoting again B ε = B 1ε B 2ε B 3ε B 4ε
, the scheme (6.101) rewrites as

1 ε R ε f n+1 j (V) f n+1 j (-V) = f n j (V) f n j (-V)) + ∆t ∆x V B 1ε j-1 2 f n j-1 (V) + B 2ε j-1 2 f n j (-V) B 3ε j+ 1 2 f n j (V) + B 4ε j+ 1 2 f n j+1 (-V) , (6.102)
where

R ε = εI 2K + ∆t ∆x V I K ζ 0 + γ + -I K ζ 0 + γ + -I K I K .
After inversion of the matrix R ε , by construction, the scheme (6.101) satisfies the well-balanced property and is asymptotic preserving. We also define:

σ 0 = K k=1 ω k e -v 2 k /2κ , σ 2 = K k=1 ω k v 2 k e -v 2 k /2κ . (6.103)
We can now state our main result for the Vlasov-Fokker-Planck equation: Theorem 6.5 The scheme (6.102) is a well-balanced approximation for the Vlasov-Fokker-Planck system (6.74). It is uniformly accurate (AP) with respect to ε. More precisely, if we assume that (6.84), (6.85) and (6.86) hold, and that σ 2 = κσ 0 in (6.103). Then, as ε → 0, the macroscopic density ρ n j := K k=-K ω k f n j (v k ) satisfies the Sharfetter-Gummel discretization (6.76).

Remark 6.6 By integration by parts, we notice that, based on (6.103),

+∞ 0 v 2 e -v 2 /(2κ) dv = κ +∞ 0 e -v 2 /(2κ) dv.
Therefore, assuming σ 2 = κσ 0 in Theorem 6.5 boils down to assume that the latter equality is also true at the discrete level.

Proof. As above, we have R ε = R 0 + εI 2K , where

R 0 := ∆t ∆x V I K ζ 0 + γ + -I K ζ 0 + γ + -I K I K .
Assuming that f admits a Hilbert expansion f = f 0 + εf 1 + o(ε), we get by identifying the terms in power of ε in (6.102),

R 0 {f 0 } n+1 j (V) {f 0 } n+1 j (-V) = 0, (6.104) and R 0 {f 1 } n+1 j (V) {f 1 } n+1 j (-V) = {f 0 } n j (V) -{f 0 } n+1 j (V) {f 0 } n j (-V) -{f 0 } n+1 j (-V) (6.105) + ∆t ∆x V B 10 j-1 2 {f 0 } n j-1 (V) + B 20 j-1 2 {f 0 } n j (-V) B 30 j+ 1 2 {f 0 } n j (V) + B 40 j+ 1 2 {f 0 } n j+1 (-V)
.

Under assumptions (6.84), (6.85), and (6.86), we may apply Lemma A.1 in Appendix. Hence

• Ker(R 0 ) = span exp(-V 2 2κ ) = span ψ 0 0 (V) , • Im(R 0 ) = Z = (Z 1 Z 2 ) ⊤ , Z i ∈ R K such that K k=1 ω k (Z 1k + Z 2k ) = 0 .
Then, equation (6.104) implies that f 0 is an element of Ker(R 0 ):

{f 0 } n+1 j (±V) = ρ n+1 j 2σ 0 exp(- V 2 2κ
).

Injecting this expression into (6.105), we deduce,

R 0 {f 1 } n+1 j (V) {f 1 } n+1 j (-V) = 1 2σ 0 (ρ n j -ρ n+1 j ) exp(-V 2 2κ ) (ρ n j -ρ n+1 j ) exp(-V 2 2κ ) (6.106) + ∆t 2σ 0 ∆x V   B 10 j-1 2 exp(-V 2 2κ )ρ n j-1 + B 20 j-1 2 exp(-V 2 2κ )ρ n j B 30
j+ 1 2 exp(-V 2 2κ )ρ n j + B 40 j+ 1 2 exp(-V 2 2κ )ρ n j+1   .
This equation admits a solution iff the right hand side belongs to Im(R 0 ). Applying again Lemma A.1, we deduce the solvability condition:

ρ n j -ρ n+1 j + ∆t 2σ 0 ∆x K k=1 ω k v k (B 10 j-1 2 exp(- V 2 2κ )) k ρ n j-1 + (B 20 j-1 2 exp(- V 2 2κ )) k ρ n j +(B 30 j+ 1 2 exp(- V 2 2κ )) k ρ n j + (B 40 j+ 1 2 exp(- V 2 2κ
)) k ρ n j+1 = 0.

Using (see (6.94) and (6.96))

γ ± exp(- V 2 2κ ) = 0 R K-1 , β ⊤ 1 R K = 1,
we deduce from the expression (6.100) of B 0 ,

K k=1 ω k v k (B 10 j-1 2 exp(- V 2 2κ )) k = 1 1 -e -E j-1 2 ∆x/κ 2E j-1 2 κ K k=1 ω k v 2 k exp(- v 2 k 2κ ); K k=1 ω k v k (B 20
j+ 1 2 exp(- V 2 2κ )) k = - e -E j-1 2 ∆x/κ 1 -e -E j-1 2 ∆x/κ 2E j-1 2 κ K k=1 ω k v 2 k exp(- v 2 k 2κ ); K k=1 ω k v k (B 30
j+ 1 2 exp(- V 2 2κ )) k = - 1 1 -e -E j+ 1 2 ∆x/κ 2E j+ 1 2 κ K k=1 ω k v 2 k exp(- v 2 k 2κ ); K k=1 ω k v k (B 40
j+ 1 2 exp(- V 2 2κ )) k = e -E j+ 1 2 ∆x/κ 1 -e -E j+ 1 2 ∆x/κ 2E j-1 2 κ K k=1 ω k v 2 k exp(- v 2 k 2κ
).

Thus, we may rewrite the limiting scheme as

0 = ρ n j -ρ n+1 j + σ 2 ∆t σ 0 ∆x   E j-1 2 κ ρ n j-1 -e -E j-1 2 ∆x/κ ρ n j 1 -e -E j-1 2 ∆x/κ - E j+ 1 2 κ ρ n j -e -E j+ 1 2 ∆x/κ ρ n j+1 1 -e -E j+ 1 2 ∆x/κ   .
By assumption σ 2 = κσ 0 , we recognize the Sharfetter-Gummel scheme (6.76).

Conclusion and outlook

In this paper, we proved (at the price of heavy technicalities) that most the discretizations advocated in [24, Part II], for 1 + 1 kinetic models endowed with a non-trivial diffusive limit, converge toward Il'in/Scharfetter-Gummel's "exponential-fitting" scheme. Such a property implies that IMEX discretizations (5.69) and (6.102), which rely on S-matrices, should be "uniformly accurate" (in the sense of, e.g., [START_REF] Berger | An analysis of a uniformly accurate difference method for a singular perturbation problem[END_REF][START_REF] Gartland | On the uniform convergence of the Scharfetter-Gummel discretization in one dimension[END_REF][START_REF] Roos | Robust numerical methods for singularly perturbed differential equations. Convection-diffusion-reaction and flow problems; 2nd edition[END_REF]) with respect to the Peclet number,

Pe = |E j-1 2 | for (5.59), Pe = E j-1 2 κ
for (6.76).

In order to ensure strong stability properties for the overall discretization, the S-matrices should be endowed with left/right-stochastic properties, as recalled in Lemma 2.5: in particular, the requirement (2.17) is important for mass-preservation. This is closely related to "matrix balancing techniques" and Sinkhorn's algorithm, [START_REF] Hershikowitz | Classifications of nonnegative matrices using diagonal equivalence[END_REF][START_REF] Sinkhorn | A relationship between arbitrary positive matrices and doubly stochastic matrices[END_REF]. Especially, results in [START_REF] Hershikowitz | Classifications of nonnegative matrices using diagonal equivalence[END_REF] ensure that one can often adjust a slightly noisy S-matrix to recover a closely related one satisfying both (2.17) and the "right-stochastic" criterion. Besides, computations in Section 6 indicate the importance of having steady solutions R + * ∋ v → Ψ ±n (x, v) being a T -system for any x > 0 and n < N ∈ N. Such an issue appears to be specific to kinetic models like Fokker-Planck equations, involving a local differential operator as the collision mechanism. Indeed, the corresponding stationary boundary-value problem is usually solved by means of a modulation of Sturm-Liouville eigenfunctions, for which the so-called "Haar property" (2.6) is not obvious; for instance it may lead to "generalized polynomials", like e.g. exponential ones, for which the equivalent condition (2.7) is not satisfied in general, as recalled in our Appendix 3.1

We now emphasize two possible applications of these schemes, one partly studied, the other being essentially an outlook extending the present work:

• Some 1 + 1 kinetic models of chemotaxis dynamics exhibit (in sharp contrast with diffusive Keller-Segel approximations) bi-stability phenomena of traveling waves, within certain ranges of parameters. Practical computational simulations were achieved in [START_REF] Calvez | Traveling chemotactic aggregates at mesoscopic scale and bi-stability[END_REF] thanks to numerical algorithms based on both S-matrices and L-splines discretizations of diffusive equations [START_REF] Gosse | Viscous equations treated with L-splines and Steklov-Poincaré operator in two dimensions[END_REF].

• In [START_REF] Cercignani | Mathematical methods in kinetic theory[END_REF]Chap. 6] and [12, Chap. 5] Cercignani explains how the formalism of "Case's elementary solutions" can be extended to linearized BGK models of the Boltzmann equation in 1 + 1 dimensions (see also the paper [START_REF] Kriese | Elementary solutions of coupled model equations in the kinetic theory of gases[END_REF] and numerical computations in [START_REF] Gosse | Well-balanced schemes using elementary solutions for linear models of the Boltzmann equation in one space dimension[END_REF]). In particular, in [12, pp. 106-108], a (formal) passage from kinetic heat transfer system toward linearized Navier-Stokes-Fourier equations. The discretization of the heat transfer system was presented in e.g. [START_REF] Gosse | Computing Qualitatively Correct Approximations of Balance Laws[END_REF]Chap. 14], so that some techniques developed in this paper may shed light onto the corresponding macroscopic behavior.

(ii) From (i), since the values of both K ∈ N and λ's are arbitrary, the set (φ λ (v)) λ≥0 clearly constitutes a Markov system.

(iii) The point (i) implies that the matrix 1 R K φ λ 1 (V) • • • φ λ K-1 (V) is invertible.

For our next property, we consider two sets of positive numbers 0 < λ 1 < . . . < λ K-1 and 0 < µ 1 < . . . < µ K-1 with corresponding Case's eigenfunctions (φ λ ) λ and (φ µ ) µ . We denote γ 1 , respectively γ 2 , the corresponding matrices defined in Proposition A.1 (iii) for the set (λ i ) i , respectively (µ i ). We introduce

ζ 1 = φ λ 1 (V) -φ λ 1 (-V), . . . , φ λ K-1 (V) -φ λ K-1 (-V) ζ 2 = φ µ 1 (V) -φ µ 1 (-V), . . . , φ µ K-1 (V) -φ µ K-1 (-V)
Then, let us denote

H = I K ζ 2 γ 2 -I K ζ 1 γ 1 -I K I K .
The following Lemma shed light onto the kernel and the range of H:

Lemma A.2 With the above notations, let us assume moreover that

∀ i, K k=1 ω k v k (φ λ i (v k ) -φ λ i (-v k )) = K k=1 ω k v k (φ µ i (v k ) -φ µ i (-v k )) = 0.
Then, the matrix H is such that:

• Ker(H) = Span(1 R 2K ),

•

Im(H) = Z = (Z 1 Z 2 ) ⊤ , Z i ∈ R K such that K k=1 ω k (Z 1k + Z 2k ) = 0 .
Proof.

• Pick Y = (Y 1 Y 2 ) ⊤ ∈ Ker(H), then Y 1 -Y 2 = ζ 1 γ 1 Y 1 = -ζ 2 γ 2 Y 2 .
Since, from Proposition A.1, the families {1 R K , φ λ 1 (V), . . . , φ λ K-1 (V)} and {1 R K , φ µ 1 (V), . . . , φ µ K-1 (V)} are basis of R K , we may write

Y 1 = a 0 + K-1 ℓ=1 a ℓ φ λ ℓ (V), Y 2 = b 0 + K-1 ℓ=1 b ℓ φ µ ℓ (V).
By definition of ζ i and γ i , i = 1, 2, we have (1λ i v)

ζ 1 γ 1 Y 1 = K-1 ℓ=1 a ℓ (φ λ ℓ (V) -φ λ ℓ (-V)), ζ 2 γ 2 Y 2 = K-1 ℓ=1 b ℓ (φ µ ℓ (V) -φ µ ℓ (-V)).
K-1 j=1

(1 + µ j v)

+ K-1 ℓ=1 a ℓ K-1 i=1,i =ℓ (1 -λ i v) K-1 j=1 (1 + µ j v) - K-1 ℓ=1 b ℓ K-1 i=1 (1 -λ i v) K-1 j=1,j =ℓ (1 + µ j v).
This is a polynomial of degree 2(K -1) which admits the 2K roots, ±v 1 , . . . , ±v K (from above equalities). So it is the null polynomial. Picking the values v = 1/λ ℓ and v = -1/µ ℓ , ℓ = 1, . . . , K -1, we deduce that a 0 = b 0 , a ℓ = 0, and b ℓ = 0, for ℓ = 1, . . . , K -1. Therefore,

Y 1 = Y 2 = a 0 1 R K .
• Consider an element in the range of H, Z = HY , with

Z = (Z 1 Z 2 ) ⊤ , Y = (Y 1 Y 2 ) ⊤ , Z i ∈ R K , Y i ∈ R K , i = 1, 2. Then, K k=1 ω k (Z 1k + Z 2k ) = K k=1 ω k v k K ℓ=1 (ζ 1 γ 1 ) kℓ Y 1ℓ + (ζ 2 γ 2 ) kℓ Y 2ℓ .
Applying our assumption, we get ∀ℓ,

K k=1 ω k v k (ζ 1 γ 1 ) kℓ = 0, K k=1 ω k v k (ζ 2 γ
2 ) kℓ = 0, so, for any Z = (Z 1 Z 2 ) ⊤ ∈ Im(H), we have K k=1 ω k (Z 1k + Z 2k ) = 0. The dimension of Ker(H) is 1, so, thanks to the rank-nullity Theorem, equalities are as claimed in Lemma A.2.

Properties of eigenelements of VFP

This appendix is devoted to the proof of an analogue of Lemma A.2 for the VFP case under assumptions on the set of discrete velocities. We first define the useful notations. Let ψ 0 ℓ , ℓ = 0, . . . , K -1, be defined as in (6.83). Let us assume that assumptions (6.84), (6.85) and (6.86) on the velocity quadrature hold. Therefore, there exists β ∈ R K and γ ∈ M K-1×K (R) such that

β ⊤ γ = ψ 0 0 (V) ψ 0 1 (V) • • • ψ 0 K-1 (V) -1
.

We introduce ζ ℓ := ψ 0 ℓ (V)ψ 0 ℓ (-V), and ζ := ζ 1 . . . ζ K-1 ∈ M K×K-1 (R). Then we denote the matrix

H = I K ζγ -I K ζγ -I K I K .
Lemma A.1 With the above notations, if we assume that (6.84), (6.85), and (6.86) hold. Then, we have

• Ker(H) = span exp(-V 2 2κ ) = span ψ 0 0 (V) ,

• Im(H) = Z = (Z 1 Z 2 ) ⊤ , Z i ∈ R K such that K k=1 ω k (Z 1k + Z 2k ) = 0 .

Proof. We proceed as in the proof of Lemma A.2.

• Let Y = (Y 1 Y 2 ) ⊤ ∈ Ker(H), then

Y 1 -Y 2 = ζγY 1 = -ζγY 2 .
By assumption (6.84), the family {ψ 0 0 (V), ψ 0 1 (V), . . . , ψ 0 K-1 (V)} is a basis of R K , then, we may write

Y 1 = K-1 ℓ=0 a ℓ ψ 0 ℓ (V), Y 2 = K-1 ℓ=0 b ℓ ψ 0 ℓ (V).
Simple computations using the definition of ζ and γ and recalling that ψ 0 0 (V) = exp(-V 2 2κ ), give From assumption (6.85), we deduce that a 0 = b 0 , a ℓ = 0, b ℓ = 0, for ℓ = 1, . . . , K -1. As a consequence Y 1 = Y 2 = a 0 ψ 0 0 (V). • Consider an element in the range of H, Z = HY , with

ζγY 1 = K-1 ℓ=1 a ℓ (ψ 0 ℓ (V) -ψ 0 ℓ (-V)), ζγY 2 = K-1 ℓ=1 b ℓ (ψ 0 ℓ (V) -ψ 0 ℓ (-V)).
Z = (Z 1 Z 2 ) ⊤ , Y = (Y 1 Y 2 ) ⊤ , Z i ∈ R K , Y i ∈ R K , i = 1, 2. Then, K k=1 ω k (Z 1k + Z 2k ) = K k=1 ω k v k K ℓ=1 (ζγ) kℓ Y 1ℓ + Y 2ℓ .
Applying our assumption, we get ∀ℓ, K k=1 ω k v k (ζγ) kℓ = 0, so, for any Z = (Z 1 Z 2 ) ⊤ ∈ Im(H), we have K k=1 ω k (Z 1k +Z 2k ) = 0. The dimension of Ker(H) is 1, so, thanks to the rank-nullity Theorem, rank(H) = K -1, which allows to conclude the proof.

Remark 4 . 2 = 1 2 K

 422 Under assumption (4.31) on the discrete velocity set, we may choose the measure ν k=-K ω k δ(vv k ) in this proposition. Then, relation (4.34) becomes

Lemma 6 . 4

 64 For E > 0,B 0 := lim ε→0 B ε = B 10 B 20 B 30 B 40 ,(6.100)

Thus from the equalities Y 1 = 1 ℓ=1a 1 ℓ=1a

 111 Y 2ζ 2 γ 2 Y 2 and Y 2 = Y 1ζ 1 γ 1 Y 1 , we deduce a 0b 0 + Kℓ φ λ ℓ (V)b ℓ φ µ ℓ (-V) = 0, a 0b 0 + Kℓ φ λ ℓ (-V)b ℓ φ µ ℓ (V) = 0.We now proceed as in the proof of Proposition A.1 by introducing the polynomial v → Q(v) := (a 0b 0 )

Thus the equalities Y 1 = Y 2 -a

 12 ζγY 2 and Y 2 = Y 1 -ζγY 1 imply (a 0b 0 ) exp(ℓ ψ 0 ℓ (-V)b ℓ ψ 0 ℓ (V) = 0.

  Let n ∈ N and F n = (f 1 , f 2 , ..., f n ) be a family of functions, continuous on an interval I ⊂ R: it is endowed with the Haar property if, for any strictly increasing family X

	We first recall basic notions from standard (one-dimensional) approximation theory, following mostly
	[15, Chapter 3].
	Definition 2.1

2.1 Haar property, Chebyshev T -systems and Markov systems

  .[START_REF] Scharfetter | Large signal analysis of a silicon read diode oscillator[END_REF] We are now in position to state our main result for the radiative transfer equation, whose proof is postponed to the next subsection. The scheme (4.49) for the radiative transfer equation (4.29) is well-balanced and uniformly accurate (AP) with respect to ε. Moreover, assuming (4.31) holds, when ε → 0, the macroscopic density ρ n j := K k=-K ω k f n j (v k ) solves the centered finite difference scheme for the heat equation (4.30).As notice in §2.3, we may rewrite (4.49) in IMEX form (see(2.14)),

	Theorem 4.6

Acknowledgment

We gladly thank Prof. Christian Krattenthaller (Vienna) for his kind help in the study of the Haar property satisfied by exponential monomials.

This work is supported by French/Italian PICS project MathCell (CNRS/CNR). NV acknowledges partial support from french "ANR Blanche" project Kibord: ANR-13-BS01-0004.

Appendix 1 Some properties on Case's eigenelements

In this Appendix, we establish some useful properties on the eigenfunctions defined in Proposition 4.1. First we show that the set of Case's eigenfunctions is endowed with the Haar property (see Definition 2.1).

Proposition A.1 Denote φ λ (v) = 1 1-λ v for λ ≥ 0, then the following properties hold:

(ii) The set (φ λ ) λ≥0 is a Markov system on R + * in the sense of Def. 2.2.

(iii) There exists β ∈ R K and γ ∈ M K-1×K (R) such that

.

Proof. These properties are shown by studying convenient polynomials.

(i) As the family contains K vectors, it suffices to show its linear independence: Assume

Thus, for any k ∈ {1, . . . , K}, the polynomial

has degree K -1, but K roots {v 1 , . . . , v K }, so that P is a null polynomial, i.e.

∀v ∈ R, P (v) = 0.

Identifying the term of higher degree, we deduce that a 0 = 0. Then, taking v = 1/λ i , i ∈ {1, . . . , K -1}, we obtain a i = 0, for all i ∈ {1, . . . , K -1}.

3 Some properties of exponential polynomials

Elementary proof of the Pólya-Szegö estimate

Hereafter, following [37, page 10], we establish by induction a simple bound on the number of real roots of an exponential polynomial; for various extensions, see [START_REF] Voorhoeve | On the Oscillation of Exponential Polynomials[END_REF][START_REF] Wielonsky | A Rolle's theorem for real exponential polynomials in the complex domain[END_REF] ∀n ∈ N,

We intend to show that, for any n ∈ N, f n admits at most N n -1 roots, where

We use an induction on n:

• for n = 1, the exponential polynomial reads f 1 (x) = P 0 (x) exp(µ 0 x) so it admits at most k 0 = N 1 -1 roots.

• Assume the property (A.1) holds for f n , so that it admits at most N n -1 real roots. Let M be the number ot real roots of f n+1 , and define

By the classical Rolle's theorem for smooth functions, its (1

admits at least M -(1 + k n ) roots. But since g n+1 is an exponential polynomial to which (A.1) applies, it comes that

Haar property for exponential monomials

Although the former estimate suggests that exponential polynomials is not a Chebyshev T -system, non-negative exponential monomials do satisfy the Haar property on (0, +∞):

Theorem A.1 (Krattenthaller, [START_REF] Krattenthaler | [END_REF]) Let (x 0 , x 1 , . . . , x n-1 ) ∈ R n + , (y 0 , y 1 , . . . , y n-1 ) ∈ R n + be non-negative with y 0 < y 1 < • • • < y n-1 . Moreover, let (z 0 , z 1 , . . . , z n-1 ) ∈ N n be non-negative integers with z 0 < z 1 < • • • < z n-1 . The generalized Vandermonde determinant, det 0≤i,j<n e y j x i x z j i (A.2) vanishes if and only if two of the x i 's are equal to each other.

The proof of this result relies on an expansion of the exponential and the use of Schur functions [START_REF] Krattenthaler | Watermelon configurations with wall interaction: exact and asymptotic results[END_REF][START_REF] Lindström | On the vector representations of induced matroids[END_REF][START_REF] Macdonald | Symmetric Functions and Hall Polynomials[END_REF].


Philippe Dessus
email: philippe.dessus@upmf-grenoble.fr

Stefan Trausan-Matu
email: stefan.trausan@cs.pub.ro

Sonia Mandin
email: sonia.mandin@msh-alpes.fr

Traian Rebedea
email: trebedea@gmail.com

Virginie Zampa
email: virginie.zampa@gmail.com

Mihai Dascalu
email: mikedascalu@yahoo.com

Emmanuelle Villiot-Leclercq
email: emmanuelle.villiotleclercq@gmail.com

In Workshop "Analysing the quality of collaboration in task-oriented computer-mediated interactions"

We describe the main tasks students usually complete when working in an elearning platform, across five mean features that have to be taken into account in research efforts (writing-based activities, individual/collective level, knowledge/ pedagogy orientation, feedback, multiple stakeholders account). Ways to analyse and assist these tasks by (semi)-automatic assessments using NLP techniques is eventually discussed.

Introduction

Current e-learning platforms allow rich collaborative learning activities that are now very well detailed and documented [START_REF] Dillenbourg | Over-scripting CSCL: The risks of blending collaborative learning with instructional design[END_REF][START_REF] Kollar | Collaboration scripts -A conceptual analysis[END_REF]. However, the ways to record, study and analyse these activities yield methodological issues often debated in the literature [START_REF] Strijbos | Methodological challenges for collaborative learning research[END_REF] and theoretical frameworks to tackle these issues are lacking. The learning activities engaged in collaborative e-learning contexts share some specificities. First, they are based on writing. Second, their manifestations are both at individual and collective level. Third, their aim is twofold: at covering (learning) a given know-ledge domain but also at leading a pedagogy-related activity. Fourth, they require analysis in order to provide an adequate feedback. Fifth and last, the stakeholders to be considered are not only the learners and the teachers, but also the researchers studying the activity. Taking into account all these specificities requires to devise ad hoc methodologies and overcoming research challenges. [START_REF] Strijbos | Methodological challenges for collaborative learning research[END_REF] listed five main methodological challenges close to those we shall point out in this paper, the goal of which is to present a comprehensive framework drawn from Bakhtin's work and a set of NLP-based tools that can help analyze learners' tasks according to these five points. The following sections shed light on each of them.

2 Five E-learning Specificities and Task Features

The tasks every learner performs in an e-learning platform share five features.

Writing for learning. Every learner engaged in individual and collaborative learning in a virtual platform performs a set of writing-based activities (e.g., abstract writing, note taking, chatting, writing in forums), which are both evidences for, and products of, learning [START_REF] Emig | Writing as a mode of learning[END_REF]. We can integrate the different writing-based learning activities in a comprehensive framework, based on Bakhtin's dialogism theory [START_REF] Bakhtin | The dialogic imagination: Four essays[END_REF]. As [START_REF] Koschmann | Toward a dialogic theory of learning: Bakhtin's contribution to understanding learning in settings of collaboration[END_REF] put it, quoting Bakhtin: "[…] the voices of others become woven into what we say, write and think". We thus can take into account all these activities within a unique framework: everything-written, read or spoken-has a dialogic nature, which is expressed through writing and relates to learning.

Multilevel Tasks: from Individual to Collective. Tasks carried out by students are often separated in two independent ones, individual and collective. As [START_REF] Stahl | Group cognition. Computer support for building collaborative knowledge[END_REF] puts it, learners engaged in a collaborative task in a e-learning platform have to cope with two recursive and interrelated main tasks: first, they are involved in an individual knowledge-building process; second, they are publicly engaged in a process of collaborative discussions about the notions at hand in the first loop. Bakhtin's ideas of dialogism and inter-animation suit with these intertwined and multiple tasks (inner dialogs and debates).

Two aims: knowledge and pedagogy. The multiple tasks in which students are engaged in e-learning do not share the same goals. The complexity of any learning situation is partly due to the fact that two different and often conflicting aims interact with each other [START_REF] Shulman | Those who understand: Knowledge growth in teaching[END_REF]): learning a knowledge domain and in parallel being confronted with pedagogy-driven activities. On the side of knowledge, learners are given information they process in order to acquire knowledge. On the side of pedagogy, learners' behavior is directed as 'moves' within the classroom environment and pedagogical methods can be inferred from these moves.

Feedback delivery.

In an e-learning context, students spend lot of time waiting for feedback from teachers or tutors about their writing, whatever are the goals and levels pursued. They encounter some problems: they stagnate themselves in the writing process; the limited feedback opportunities do not stimulate explorative approaches ("what if-trials"), but force them to hand in mainly completed versions; during writing, it is difficult to self-assess ongoing work and understanding. Teachers have a limited overview of the learners' processes, and assessments of students' understanding or collaboration are difficult and timeconsuming. Feedback is thus necessary in e-learning contexts and can partly be automated by computer-based procedures.

Accounting for stakeholders' viewpoint. E-learning contexts are populated by numerous stakeholders (students, tutors, teachers, researchers) whose tasks may differ, overlap or be contradictory to each other. These tasks can also strongly interfere with the kind of tool used for analysing a given learning situation. Since most of the tools aiming at analysing collaborative software are devised for research purposes, they are more difficult to be used by other stakeholders.

NLP-Based Tools

Web-based services using NLP techniques can take into account the five features of e-learning situations presented above:

(1) detection of relations between utterances can be processed to reveal the voices engaged in writing or dialog; (2) account for both the individual and collective level of knowledge acquisition;

(3) sensitivity to both knowledge (cognitive models) and 'moves' (dynamic situations) [START_REF] Dessus | L'enseignement en tant que supervision d'un environnement dynamique[END_REF][START_REF] Wolfe | Learning from text: Matching readers and texts by Latent Semantic Analysis[END_REF]); (4) possibility to deliver just-in-time feedback allowing self-paced learning;

(5) deliver generic feedback to account for all the stakeholders' categories. Let us now present two instances of web-based services designed from this viewpoint. Pensum supports learners at an individual level in the automatic assessment of their essays (summaries, syntheses). Pensum analyses how well learners understand course texts through their textual productions. It provides different kinds of feedback (see Figure 1) all based on LSA (Latent Semantic Analysis, [START_REF] Landauer | A solution to Plato's problem: the Latent Semantic Analysis theory of acquisition, induction and representation of knowledge[END_REF] on two important features influencing writing quality: topic coverage (semantic links between sentences source texts and synthesis) and inter-sentence coherence.

PolyCAFe (Chat & forum Analysis and Feedback system, Trausan-Matu & Rebedea 2010) functions at a collective level using an NLP pipe (stemming, POS tagging, chunking, etc.), advanced pattern matching, social network analysis and LSA for detecting discussion topics, threads and inter-animation in chat logs.

Feedback (textual and graphical) is generated emphasizing collaboration degree, discussed topics and evaluation of the participant's contribution (see Figure 2).

Conclusion and Future Work

We presented a theoretical framework in which five main learning tasks occurring in e-learning collaborative platforms, as well as two web-services fitting with these tasks. Our services (1) are focused on writing activities; (2) are both on individual and collective levels (3) can embed pedagogical facets through the use of widgets; (3) propose high-level and automated feedback; (5) can be used by different stakeholders. Further work is planned to cross the results of these tools to uncover patterns of efficient individual or collaborative forms of writing.

Figure 1 .

 1 Figure 1. The different pieces of feedback delivered by Pensum.

Figure 2 .

 2 Figure 2. PolyCAFe main interface.

Acknowledgments

This work was partly supported by the 7th Framework Programme (Theme FP7-ICT-2007-1-4.1) of the European Community (LTFLL-STREP Project, Grant no 212578).