On the numerical resolution of an ion / electron collision operator of the Fokker-Planck type coupled with an electronic equation C. Buet1, S. Dellacherie2, R. Sentis1 Commissariat à l'Énergie Atomique Centres d'étude de 1Bruyères-le-Châtel et de 2Saclay France Abstract

By introducing the denition of the entropic average, we propose a new numerical scheme which solves the collision operator ion / electron coupled with an electronic equation. We show that this new scheme, which can be explicit or semi-implicit, has good theoretical properties and that the numerical results are good even on a rough grid.

Introduction We consider a one species hot homogeneous plasma and we only study the electron / ion interactions. The time scale is then the ion / electron collision time. The model consists of a Fokker-Planck equation for the ions, a macroscopic equation for the electronic temperature and the quasi-neutrality relations. We write the model in the following way : 8 > > @t f (t;v) = S (f ); > > > > > < > > > > > > > :

@t ( 2d NeTe ) = d N (Te T );

Ne = ZN; U!e = ! U < ! v f > =N

(1)

with N =< f > and 2d NT =< 12 m(v ! U )2 f > R (see [1]) where < g > g(v)d! v in the continuous case and < g > P g ( v )
v in the discrete case, d being the dimension of the velocity j j

1

21st International Symposium on Rareed Gas Dynamics space. f (t; v) is the ionic distribution, (Te ) = C ste=Te3=2 is the collision frequency and Te is the electronic temperature; m and Z are respectively the atomic mass and the atomic number. S (f ) describes the collisions between ions and electrons and is dened by h i S (f )  rv
(! v U!e )f (t; v) + Tme rv f (t; v)

in the convection-diusion form. Denoting

Mf (t; v) =

N exp (2Te =m)d=2



S (f ) can be written in the Landau form 

 m(! v U!e )2 ; 2Te 

S (f )  mTe rv
f (t;v)rv log Mf (t;(t;v)v) f



(2)

whose discretization naturally leads to a numerical version of the H theorem:   f ( t;v ) < log M (t; v) S (f ) > 0 and S (f ) = 0 () f (t;v) = Mf (t; v): f On the other hand, it is easy to verify that we have the following macroscopic conservation equations : 8 < :

@t < f >= @t N =!0; @t < ! v f >= N@t U = !0 ; @t < 12 mv2 f >= d N (Te T ):

The plan of this article is the following : In the rst section, we study the semi-discrete numerical scheme established from

the convection-diusion form of S (f ) to solve the non linear system (1), and we introduce the denition of the entropic average to evaluate the ionic distribution f at the interfaces of the mesh. This study is important because it allows us to prove that there exists a unique maxwellian equilibrium state for the semi-discrete, the explicit and the semi-implicit schemes which minimizes the Boltzmann entropy in the case of the semi-discrete and explicit schemes. In the second section, we establish the explicit numerical scheme coming from the semi-discrete numerical scheme. This explicit scheme preserves the mass, the impulsion, the energy of the plasma and the thermodynamical equilibrium. Moreover, we prove that, under a CFL criterium, this scheme preserves the positivity of the density f and of the electronic temperature Te and that the entropy decreases. In the third section, we describe the semi-implicit version of this explicit scheme. This scheme preserves the mass, the impulsion, the energy and the thermodynamical equilibrium too. In the fourth section, we show with numerical tests that this new semi-implicit scheme has good

2

Marseille, France, July 26-31, 1998 properties similar to those found from a discretization of the Landau form (see for example [2] and [3]) and that it is more precise than the classical Chang and Cooper's scheme (see [4]) when the discretization is rough. In that paper, we take the velocity dimension d equal to 1, the extension in the cartesian case to higher dimension being easy. The demonstrations of the propositions and theorem written in this paper are done in [6].

1 On the semi-discrete numerical scheme Let's dene the following semi-discrete numerical scheme 8 < @t fj = S (f; Te )j ;

(3)




@t 21 NeTe = (Ne Te NTe )

:

where v 2 Rd=1 with h i S (f; Te )j =
v (vj+1=2 Ue )fej+1=2 (vj 1=2 Ue )fej 1=2 + m
Tve 2 (aj fj+1 bj fj + cj fj 1 ) where P 8 Ne = fej+1=2
v; > > > > > > > <

(4)

j

P

vj+1=2 fej+1=2
v=Ne + ue; j > > > > > P > > : Te = m(vj+1=2 Ue )2 fej+1=2
v=Ne + et and

Ue =

j

8 > <

ue = mTNe
(fjmax f1 ) ;

> :

e

h







i

et = TNe
fjmax vjmax+1=2 Ue + f1 Ue v1=2 with aj = cj = 1 and bj = 2 except at the boundary of the mesh: see the boundary limit conditions. fej+1=2 (t) is an evaluation of f (t = tn ; v = vj+1=2 ) and is the entropic average of fj (t) and of fj+1 (t) e

dened by :

Denition The entropic average fej+1=2 of fj > 0 and of fj+1 > 0 is

dened by

fej+1=2 = log ff

j j

fj+1 if f 6= f ; j j+1 log fj+1

= fj else.

3

21st International Symposium on Rareed Gas Dynamics By continuity, we extend this denition to the case where fj+1 = 0 and / or fj = 0 by writing

fej+1=2 = 0 if fj
fj+1 = 0:

ue is a correction which allows the scheme to exactly conserve the impulsion (et is a consequence of the fact that ue 6= 0).

Boundary limit conditions We use a Robin's boundary limit condition. We then have at the boundary of the mesh (v U )f + Tm @v f = 0 which is equivalent to taking for the numerical scheme 8 a = 1 if j 6= jmax; > < j bj = 2 if j 2 f2; :::;jmax 1g ; > : cj = 1 if j 6= 1; b1 = bjmax = 1 et ajmax = c1 = 0 and fe1=2 = fejmax+1=2  0: This Robin's boundary condition coupled with the correction u~ gives the following conservation property : e

Property 1.1

8 > > > > < > > > > :

@t < fj >= @t N = 0; @t < vj fj >= N@t U = 0; 2

@t < m v2 fj >= (NTe Ne Te): j

Remark that this property would be veried by using for example the arithmetic average instead of the entropic average. Let's now dene

H (f;Te )(t) =< fj log fj > ZN 2 log Te ; 



e 2 Mf;T (t; v) = p N exp m(v2TeU ) 2Te =m   N m(v Ue 1 )2 p and M1 ( v ) = exp f;T 2Te1 2Te1 =m where Ue 1 and Te1 are dened with the initial conditions U 0 , T 0 and Te0 in the following way : 8 0 1 < < (vj U )
Mf;T ;j >= 0; e

e

e

:

1

0

< f[ 21 m(vj U 0 )2 + Z T2 ] [ T20 + Z T2 ]g
M1 f;T ;j >= 0: e

e

We then have the following proposition :

4

e

(5)

Marseille, France, July 26-31, 1998

Proposition 1.1 The system (5) has an unique solution M1f;T . e

Let's nally dene 0 f 1 = < MN1 > M1 f;T f;T ;j

e

e

0

1 and H (f 1; Te1 ) =< fj1 log fj1 > ZN 2 log Te : The main properties of the scheme are summerized in the following theorem:

Theorem 1.1 For all strictly positive initial condition, when fej+1=2 (t)

is the entropic average of fj (t) and of fj+1 (t), the semi-discrete scheme dened with (3) and (4) assures that : i) @t H (f;Te )  0; inf Te (t) > 0 and

t2[0;+1[

8t : H (f 1 ; Te1 )  H (f; Te )(t)

(Gibbs lemma): ii) The series fj (t) converges in the l1 norm towards the unique thermodynamical equilibrium fj1 and Te (t) converges towards Te1 when t goes to the innity.

2 On the explicit numerical scheme

Let's now dene the explicit numerical scheme (with again v 2 Rd=1 ): 8 1 n+1 f n ) = S (f n ; T n )j ; > j e <
t (fj i h
> 1 NeTe
n = n(N e n Ten :
1t 21 Ne Te n+1 2

(6)

N n Ten )

where S (f; Te )j , Ne Te and the boundary conditions were dened in the previous section. Let's note that this scheme veries the discrete time version of Property 1.1. We now dene hnmin
v 2 Z
1 and
tn2 = 4


tn1 = 4 m n Ten
hnmax
Mn (1 + n ) n 1 "n hM

with Mn = maxj M

n ;j f n ;Te

1

n ;j f n ;Te

and

hnmax = max k



fn

Mf

n

i

n

n



;Ten k

;

We then have the following proposition :

5

n n e

hnmin = min k

4

v n = Zmax (T =m)2

"n = NNe
Tet ;

;

k n e



fn

Mf

n



;Ten k

:

k

21st International Symposium on Rareed Gas Dynamics

Proposition 2.1 For all strictly positive initial condition, when fejn+1=2 n n

is the entropic average of fj and of fj+1 , the explicit numerical scheme dened with (6) preserves the positivity of fjn+1 and of Ten+1 and veries H (f 1 ; Te1 )  H (f n+1 ; Ten+1 )  H (f n; Ten ) when
t < min(
tn1 ;
tn2 ): (7)

3 On the semi-implicit numerical scheme

Let's dene the following semi-implicit numerical scheme : 8 1 n+1 f n ) = S (f n ; f n+1 )j ; > j <
t (fj > :

with

1
t

h


n+1 1 2 Ne Te


n 1 2 NeTe

nh

i

= n(Ne n Ten N n Ten )

(v e n en S (f n ; f n+1 )j =
v j+1=2 U )fj+1=2 (vj

1=2

(8)

Ue n )fejn 1=2

i

(9) n Ten

n +1 n +1 n +1 + m
v2 (aj fj+1 bj fj + cj fj 1 ): Note that, as in the explicit case, this scheme veries the discrete time version of Property 1.1. We now have the following proposition:

Proposition 3.1 When fejn+1=2 is the entropic average of fjn and of n

fj+1 , the semi-implicit numerical scheme dened with (8) and (9) veries : f n = f 1 and Ten = Te1 () f n+1 = f n and 8j : fjn > 0; Ten > 0:

4 Numerical results

In this section, we present numerical results concerning the semiimplicit numerical scheme. The initial conditions are such that the ionic temperature and the electronic temperature are respectively 1 Kev and 2 Kev, the initial ionic distribution being a maxwellian. We take Z = 2 : we then have T = Te1 = 35 Kev when t = +1 in the continuous case. On the gures 1 to 4, we can see that the entropic average gives very good and precise results when the mesh is not thin, which is not the case if we replace the entropic average with the arithmetic average. Moreover, we can see that the entropic average is more precise than the classical Chang and Cooper average (cf. [4]). We can verify too that this semi-implicit scheme accepts very higher time steps than the totally explicit scheme.

6

Marseille, France, July 26-31, 1998

Conclusion

In this paper, by introducing the denition of the entropic average, we presented a new numerical scheme to solve the homogeneous nonlinear ion / electron collision operator of the Fokker-Planck type coupled with the macroscopic electronic equation. We showed that the semi-discrete scheme has the same properties than the continuous model in the cartesian geometry (i.e. convergence towards a unique maxwellian thermodynamical equlibrium and existence of this equilibrium). After that, we showed that the explicit version of this numerical scheme preserves the positivity of the ionic distribution and of the electronic temperature, and that the Boltzmann's entropy decreases under a CFL criteria. Finally, we showed that the semi-implicit version of this scheme preserves the thermodynamical equilibrium when it is reached. Then, we presented numerical results which shows that this new semi-implicit numerical scheme is very precise even when the mesh is not thin, property which is very important if we want to solve the Fokker-Planck equations on a grid with not too much meshes. Finally, it is important to say that all these results can be extended to the axisymmetric case (cf. [6]), property which is important if we want to simulate the Fokker-Planck equations in the case of the Inertial Connement Fusion.

References [1] M. Casanova, O. Larroche and J.P. Matte - Kinetic Simulation of a Collisional Shock Wave in a Plasma - Physical Review Letter, Volume 67, Number 16, October 1991. [2] Yu. A. Berezin, V. N. Khudick and M. S. Pekker - Conservative Finite-Dierence Schemes for the Fokker-Planck Equation Not Violating the Law of an Increasing Entropy - Journal of Computational Physics, Volume 69, Number 1, March 1987. [3] C. Buet and St. Cordier - Numerical Analysis of Conservative and Entropy Schemes for the Fokker-Planck-Landau Equation To appear in SIAM, 1999. [4] J.S. Chang and G. Cooper - A Practical Dierence Scheme for Fokker-Planck Equations, Journal of Computational Physics, Volume 6, pp. 1-16, 1970. [5] V.A. Mousseau and D.A. Knoll - Fully Implicit Kinetic Solution of Collisional Plasmas - Journal of Computational Physics, Volume 136, pp. 308-323, 1997. [6] S. Dellacherie - Contribution à l'analyse et à la simulation numérique des équations cinétiques décrivant un plasma chaud Ph.D. Thesis of the Denis Diderot University (Paris VII), 1998.

7

21st International Symposium on Rareed Gas Dynamics

T (Kev)

T (Kev)

Mesh : Jmax=10

Mesh : Jmax=10 and 100

Te

Te Mesh : Jmax=100

Arithmetic average

T T

Entropic average

Mesh : Jmax=10 t (ns)

t (ns)

Figure 2: T (t) and Te (t)

Figure 1: T (t) and Te (t) T (Kev)

Distribution f (t=+∞) ++++ + + + +

Mesh : Jmax=10 and 100

Te

+ O ×+

Entropic average (x)

+ +

and

+ +

Chang and Cooper average

T

+

Arithmetic average (o)

Maxwellian distribution (+)

+ O ×+ + + + + +

+

+

+

+

+

Jmax=10

+

+

+

+

+

+ × O + + + + + + + + O +++++++++×+++ × 0 ×O+++++++++++++ O

(x)

(o)

+ × O + + + + + ++ +++ ×++++++++++++++ O +++++++++ × × O O

Velocity (cm/s)

t (ns)

Figure 4: Distribution f (t = +1)

Figure 3: T (t) and Te (t)

8