Well-Posedness of a Diffusive Gyrokinetic Model

We study a finite Larmor radius model used to describe the ions distributions in the core of a toka-mak plasma, that consist in a gyro-kinetic transport equation, coupled with an electro-neutrality equation. Since the last equation do not provide enough regularity on the electric potential, we introduce a simple linear collision operator adapted to the finite Larmor radius approximation. Next we study the two-dimensional dynamics in the direction perpendicular to the magnetic field and prove thanks to the smoothing effects of the collisions and of the gyro-average the global existence of solutions, as well as short time uniqueness and stability.

The model studied in that article describes the density of ions in the core of a tokamak plasma. In such highly magnetized plasma, the charged particles have a very fast motion of gyration around the magnetic lines, called the Larmor gyration. A good approximation is then to consider that the particles are uniformly distributed on gyro-circle, parametrized by their gyro-center, and Larmor radius r L (that is proportionnal to the speed of rotation u, and in our article, we will forget the physical constant and write r L = u). The models obtained in that new variables are kinetic in the direction parallel to the magnetic field lines, and fluids (precisely a superposition of fluid models) in the perpendicular direction. For rigorous derivation of such models and more complete discussion on its validity, we refer to FreSon [3] and our previous work GheHauNou09 [ 4], in which the derivation is perform from a Vlasov equation in the limit of large magnetic field.
Such gyro-kinetic models are usually closed by an electro-neutrality equation, that as usual provide very few regularity for the eletric field, so that the well-posedness of gyro-kinetic models is, at least at our knowledge unknown. In this article, we add a "gyro-averaged" collision operator to the model and study the dynamics in the directions perpendicular to the field only.
Let us now describe our precise model. The ion distribution function f (t, x, u) in gyro-coordinates depends on the time t, the gyro-center position x ∈ T 2 and the velocity of the fast Larmor rotation u ∈ R + (which is also proportional to the Larmor radius). The electric potential Φ depends only on (t, x). They satisfy the following system of equation on Ω = T 2 × R + ∂f ∂t where β and ν are two positive constant, ρ is the density in physical space, T is the ion temperature, We also used the notation b ⊥ = (−b 2 , b 1 ), for any vector b = (b 1 , b 2 ) of R 2 .
That model without the Fokker-Planck operator (ν = β = 0) was studied in our previous work GheHauNou09 [4] -to which we refer for an heuristical derivation of the electro-neutrality equation ( eq:elecneutr

1.2) -and is used by physicists for simulation, by instance in Gysela code
Gysela06 [5]. Here we just mention that ( eq:elecneutr 1.2) is obtained in a close to equilibrium setting, with an adiabatic hypothesis on the distribution of the electrons n e = n 0 e− eΦ Te ≈ n 0 1 + eΦ Te , and an hypothesis of adiabatic response of the ions on the gyro-circle which gives rise to the Φ * H T term. As usual in quasi-neutral equation, we have no good a priori estimates on the regularity of E = −∇Φ.
Remark that even if that equation ( eq:gyroFP2D 1.1) is derived from a Vlasov model (A rigourous derivation of a more general 3D model is performed for fixed field E in section sect:dif 2), it is of "fluid" nature. In fact there is no transport in the remaining of the velocity variable u, and the position of the gyro-center is transported by the eletric drift (J 0 u E) ⊥ . So that the equation is similar to the 2D Navier-Stokes equation written in vorticity. More precisely, we have a family of fluid model depending on a parameter u, which are coupled thanks to diffusion in the u variable, and by the closure used for E described below.
Moreover, we will prove in the following that thanks to gyro-average J 0 u , equation has the same regularity than the NS2D equation in vorticity. In fact, the force field J 0 u ∇ x Φ belongs naturally to H 1 if f ∈ L 2 with some weight. That is why we obtain the same result that are known about the NS2D equation : global existence and short time uniqueness and stability. However, our model present an additional difficulties which is the lose of regularity for small u. In fact, for small value of u the H 1 bound (in x only) of J 0 u ∇ x Φ explodes. To state our reuslts properly, we will need the following definitions and notations : • In the sequel, the letter C will design a numerical constant, that may change form line to line. Unless it is mentioned, such constants are independent of everything.
• L 2 u (Ω) = L 2 (Ω, udxdu) is the space of square integrable functions with respect to the measure udxdu. • We shall use various norm on T 2 or on Ω. To avoid confusion, we will use the following convention. All the norm performed on the whole Ω will have their weight with respect to u as additional indice. By instance · 2πu , · H 1 2πu(1+u 2 ) . All the norms without any indices are norm on T 2 only.
• For any weight function k : R + → R + , the norm · 2,m is defined for any function f on Ω by • The most usefull weights will be m(u) = 2πu(1 + u 2 ) andm(u) = 1 + u 2 .
• We change a litlle bit the duality used to define distributions in the following definition Definition 1.1. Using distributions with the weight u means that duality is performed as This definition may seem a little artificial because the simple definition of derivative with respect to u, is not valid. Instead, However, this weight respects the underlying physics (u is in fact the 1D norm of a 2D velocity variable) and has many advantages. For instance the operator (1/u)∂ u (u∂ u ) is self-adjoint with this weight.
Our precise result are the following. We prove global existence under the hypothesis f i 2,m < +∞. :existence Then there exists at least one weak solution f ∈ L ∞ (R + , L 2 u (Ω)) ∩ L 2 (R + , H 1 u (Ω)) to ( eq:gyroFP2D 1.1)-( eq:elecneutr 1.2) with initial condition f i , which also satisfies for any t > 0 Then lim This local result has some more consequence when relating it to the bound on ν In the next section the diffusive operator of ( eq:gyroFP2D 1.1) is rigourously derived from a linear Vlasov-Fokker-Planck equation in the limit of large magnetic field. In the third section, some useful lemmas are established, proving regularizing properties of the gyro-average, global preservation of some weighted norm of f , the short time preservation of the u(1 + u 2 )-moment of ∇ x f by the system (  In that section, we rigorously justify the form of the Fokker-Planck appearing in the right-hand side of ( eq:gyroFP2D 1.1). The usual collision operator for plasmas is the nonlinear Landau operator originally introduced by Landau Landau [7]. Because of its complexity, simplified collision operators have been introduced. An important physical litterature exists on the subject, also in the gyro-kinetic case (See

Brizard04
[2] and the references therein). In this paper we choose the simplest possible operator possible, namely a linear Fokker-Planck operator. The reasons of this choice are : -Its simplicity will allow to focus on the other difficulties of the model, -The fact that physicists studying gyro-kinetic models for the core of the plasma mainly assume that the dynamics stays close to equilibrium, in which case a linear approximation of the collision operator is relevant.
-The aim of the paper is not a precise description of collisions. In fact, even if they exists in tokamaks, being needed to produce energy, their effect is small compared to the turbulent transport. However, we are interested by their regularizing effect, since the electro-neutrality equation ( eq:elecneutr 1.2) do not provide enough regularity to get a well-posed problem. This is a major difference to the Poisson equation setting.
We start from a simple model for a 3D plasma, i.e. a linear Vlasov-Fokker-Planck equation with an external electric field, an external uniform magnetic field and linear collision and drift terms, and obtain in the limit of large magnetic field a 3D (in position) equation analog to ( eq:gyroFP2D 1.1). In particular, we show that a usual linear Fokker-Planck term on the speed variables turns into an equation with diffusion terms both in space and Larmor radius variables in the limit.
Precisely, for any small parameter ǫ > 0 we study the distribution f ǫ (t, x, v) of ions submitted to an exterior electric field E(t, x) (independent of ǫ) and an uniform magnetic field B ǫ = (1/ǫ, 0, 0). We also model collisions (with similar particles and the others species) by a simple linear Fokker-Planck operator. To avoid any problem with possible boundary collisions, which are really hard to take into account in gyro-kinetic theory, we assume that (x, v) ∈ T 3 × R 3 , where T 3 is the 3D torus. When the scale length of all the parameters are well chosen (in particular the length scale in the direction perpendicular to the magnetic field should be chosen of order ǫ times the length scale in the parallel direction, we refer to our previous work GheHauNou09 [4] for more details on the scaling), the Vlasov equation f ǫ satisfies is where β, ν are two positive parameters, the subscript (resp. ⊥) denotes the projection on the direction parallel (resp. on the plane perpendicular) to B, and the superscript ⊥ denotes the projection on the plane perpendicular to B composed with the rotation of angle π/2. In others The next results require the additional notation, which is a gyro-average performed in phase space, that will be used as an initial layer to adapt the initial condition to the fast Larmor gyration. thm:FPgyro and f ǫ be a family of solutions to equation ( eq:vlares admits a subsequence that converges in the sense of distributions towards a functionf depending only on (t, x g , u = |v|, v || ) and solution to But in the gyro-coordinates this fast rotation is simply described by a rotation in v, The final diffusion appears in all dimensions except the x g || one. It does not mean that there is no regularization in that direction. Indeed, the models have diffusion in v || , which after some time regularize in the x g || direction. This mechanism is well known for the Fokker-Planck equation (see for instance bouchut [1]). However, we are not able to prove this phenomena in the non-linear setting because the electric field of the model lacks regularity. This is the reason why we will only study the 2D model.

Proof of Theorem
, a subsequence of f ǫ solutions of ( eq:vlares 2.1) without the collision term converges towards a solution of ( eq:FPlim

2.3) without the collision term.
In order to simplify the presentation, we will neglect the electric field and the parallel translation terms. To obtain the result in full generality, the only thing to do is to add the argument given in our previous work to the one given below. For the same reason, we shall also not treat the problem of initial conditions. So consider the above Vlasov Fokker-Planck equation without electric force field and parallel translation, The first step is to use the change of variables ( Since . Therefore, at least a subsequence of (f ǫ ) converges weakly to somē f ∈ L ∞ t (L 2 ). Passing to the limit in ( eq:vlaFP since all the other terms are bounded. For v = (ue iϕ , v || ) where ϕ is the gyro-phase, the previous equality means thatf is independent of the gyro-phase (in the sense of distribution and thus as a L 2 function). Equation ( eq:vlaFP 2.4) tested against a smooth function g independent of the gyro-phase writes We may also pass to the limit when ǫ tends to zero in this equation and obtain that the same equality holds for f ǫ replaced byf , considered as a function defined on Hence, for any function g independent on the gyrophase ϕ, it holds that The other terms appearing in ( eq:FP2 2.6) remain unchanged. Then, (2.7) eq:FPuphi Sincef is independent of ϕ, performing the integration in ϕ first makes the term containing ϕ vanish. So the functionf of the five variables (x g , u, v || ) satisfies It exactly means thatf satisfies the equation in the sense of distributions with weight u. It is the equation ( eq:FPlim 2.3) without parallel transport nor electric field.
If we look at solutions of this equation invariant by translation in the direction of B, we exactly get the 2D-model announced in the introduction. In fact, iff is a solution of ( eq:FPlim

2.3), then
is a solution of ( eq:gyroFP2D 1.1). Such an assumption on f is reinforced by experiments and numerical simulations, where it is observed that the distribution of ions is quite homogeneous in x .

Some useful lemmas ct:apriori
We prove here some a priori estimates useful for the proof of our theorem. In order to simplify the proof of some of the following Lemmas, we sometimes uses the following formulation of ( eq:gyroFP2D It is solution (in the sense of distribution with usual duality) of the following equation with 4D in space and velocity variables Heuristically, radial in u solution of equation ( eq:VFP4D 3.1) is a solution of ( eq:gyroFP2D 1.1). We can state for instance a precise Lemma in the case where φ is fixed and smooth.

Regularizing properties of the gyro-average.
In this section, some regularizing property of the gyro-average operato are proven. They are based on the fact thatĴ 0 ∼ k − 1 2 for large k (the precise bound are proved in It is important since the formula ( eq:elecneutr 1.2) giving the gyro-averaged potential in terms of the distribution f involves two gyro-averages, and thus a gain of one derivative for the gyro-averaged potential w.r.t. f . However, the regularizing properties of J 0 u are bad for small u, which raises difficulties. The first lemma of this section gives the regularity of the gyro-averaged potential in term of the potential Φ. The second one gives the regularity of the density ρ in terms of the distribution f . We will need the two following definitions before stating it.
For any U > 0, let F be a measurable function defined on The lemmas stating the regularity of Φ and ρ are the following. lem:phireg As a consequence, for any U > 0, Proof of Lemma using the bound Ĵ0 ∞ ≤ 1 proved in Lemma lem:boundJ A.1. It is the inequality i). For the second inequality, remark that For the third estimate of Lemma lem:phireg

3.2, remark that
and use the bound iii) of Lemma From this, we obtain The point iv) uses the previous inequalities. First remark that the norm H 1 U is also equal to Using this formulation and ii)-iii) leads to which gives the desired result and ends the proof of Lemma It follows from ( eq:calc 3.3) that for k = 0, Hence, sinceρ(0) = T 2 ρ(x) dx = 1 by mass conservation, As a consequence, the operator L −1 T maps any H s 0 , s ∈ R, into itself with norm 3.4 shows that L −1 is bounded for small T , and of order T for large T , the physical case of interest. The boundedness of the spatial domain is essential. When defined on the whole space R 2 rather than on the torus, the operator L −1 is not bounded. Its norm explodes in the low frequency range.

Proof of the Lemma
Consequently, This is the first claim of lemma lem:boundL 3.4 The function of w in the right-hand side of the previous inequality on |Ĥ(k)| is decreasing and goes from 1 4 at 0 to 0 at +∞. Consequently its minimal value are obtained for large w i.e. for small |k|, namely |k| = 1. Precisely , Since the Fourier representation of L −1 T is the multiplication by T (1 −Ĥ(k)) −1 we obtain that in any H s 0 , s ∈ R, which is the desired result.
3.3 Propagation of L 2 m and L 2 m (L 4 ) norms of f .

sec:moment
The two following lemmas will be useful in the sequel. Assume moreover that f i 2,m < +∞. Then any solution f satisfies and integrating in the x variable leads to Multiply the previous equation by k( u), where k is a smooth function on R 2 with compact support and integrate in the velocity variable u leads to By approximation, this is still true for functions k with unbounded supports. For k( u) = 1, Coming back to the 1D original quantities, it means that For k( u) =m(u), then ∆k = 4 and Therefore, g(t, u)m(u) d u ≤ g(0, u)m(u) d u + 2(2ν + β) 3.6 In order to simplify the presentation, we will first performed calculation without justifying every integration by parts and division. But once we obtain an a-priori result, we will explain the small adaptation needed to make it rigorous. First, we denote Multiply equation ( eq:gyroFP2D

1.1) by 3 sign(f )|f | 3 and integrating with respect to x leads to
(3.10) eq:alpha Now, we multiply bym(u), integrate with respect to u With the help of some integration by parts, we get that Thanks to that, the previous inequality simplify in Next we can estimate |∇ u α| in terms of γ. In fact by Hölder inequality So that the second term in the right hand side of the previous inequality is controlled up to a constant to the first one. We precisely get From which we conclude easily.
The other point not rigorously justified is the division by √ α that may be zero. However, since we have a diffusion equation, it may be proved that for t > 0, α > 0 everywhere. Or we can use a family of smooth approximation of √ ·. Or we can say that α + ǫ satisfy ( eq:alpha 3.10) with a additional term that has the good sign, so that it will satisfy ( eq:diffapprox 3.12), and we will obtain the desired inequality letting ǫ going to zero and then U going to infinity. It is well justified since the maximum principle applies there so that any solution with non-negative initial condition remains non-negative.

Short time estimate of the m-moment
of ∇ x f .

nt-grad_xf
The following lemma provide is central in the proof of the stability and uniqueness of the solution for short time. 3.1), written in 2D in u (with u = | u|), and obtain If we now multiply by t ∇ xf on the left and integrate in x, the function h defined by h(t, u) = 1 2 |∇ xf | 2 dx satisfies, (3.13) eq:gevol We may also multiply this equation bym(u) = (1 + u 2 ) and integrate it in u. We obtain after that (3.14) eq:gevol2 To go on, we need to understand a little better the matrix M (t, , and then remark that from there definition, J 0 and L −1 commute with derivation in x. So that our term may be rewritten Moreover, the H 2 m -norm of f appears in the right-hand side of ( eq:gevol 3.13). So that we may use it to control M . With a control on the H 1 norm of M , we do not get an infinite bound on M , like M ∞ ≤ C f H 2 . In that case, we will be able to use a classical tool to conclude. But this is almost true, we are in a critical case (d = 2 and p = 2) for the Sobolev imbeddings, but we still know that the square of M is exponentially integrable. Precisely, since M is of average 0, we have for all u > 0 the following Trüdinger inequality We refer to Moser [9] for a proof of that result. To estimate t ∇f M ∇f (1 + u 2 ) dxdu, we first perform the integral in x. For this, we apply the inequality . The previous inequality comes from the Legendre transform of e a , as the Young inequality, so that our application will be a log-exp analog of the Hölder inequalities. We obtain where we have used the Jensen inequality -precisely, for a function g of integral 1, g ln g ≤ ln( g 2 ) -in the last but one line, and the Sobolev imbedding from H 1 into L 4 with constant c s in the last one. The constant C may change from line to line. Using this and Jensen inequality in the previous equation, we get Remark that the fraction inside the logarithm is always greater than 1 so that the square root is well defined.
In order to get a bound on ∇ xf 2,m u , we use that Therefore, (3.16) eq:gevol3 We next use the inequality 1 + ln(x) ≤ x ǫ ǫ , valid for any ǫ ∈ (0, 1), x > 0, and get We will use the Young inequalities xy ≤ x p p + y q q , where 1 p + 1 q = 1. For the first term, with p = 2 1+ǫ and then q = 2 1−ǫ , and for the second with p = 4 3+2ǫ and then q = 4 1−2ǫ . We obtain the following bounds valid for ǫ < 1 2 . Taking into account the two constants ν and C(ǫ) we get for ǫ < 1 where C(ǫ) has only change from a numerical constant. With the notation h = ∇ xf 2 2,m , it gives That is a differential inequality with a growth faster than linear, that give a solution that may explode in a finite time. The time of explosion τ * may be bounded below by something depending only on ν, ǫ, T and h(0) = ∇ xfi 2,m . Since ǫ may be choosen arbitrairy between 0 and 1 2 , τ * depends only on ν,β, T and the initial value h 0 .
Using that in the inequality ( eq:gevol3 3.16), and turning back to the original u variable, we obtain the existence of a constant C 3 (T, ν, ∇ x f i m ) such that A bound by below for τ * . In order to simplify the next two paragraphs, we assume that β ≤ ν. A careful analysis of the term in the right hand side of the previous equation shows that if ǫ is choosen small enough so that ν ≤ C(ǫ), then for h ≤h := ν C(ǫ) 2 the dominant term is 4νh and for h ≥h the dominant term is h 6+ 8ǫ 1−2ǫ . Then, if h(0) ≥h, the explosion time is given by the equation For that equation, we get that the explosion time is larger than In the case h(0) ≤h, then for the early time, the equation may be rewritten till h(0) =h. It take a time greater than T 1 = 1 12ν ln h h(0) . And after that time, the explosion time is given by the later calculation, and due to simplification it comes C ν . Finally, we get an explosion time Best choice for ǫ. It is quite difficult to optimize that quantity in ǫ. But, as the condition on ǫ are 0 < ǫ < 1 2 and ν ≤ C(ǫ) = C cT ǫ , we can choose ǫ = min With that choice, we get else It is then clear that the value of τ * depends only on ν, h 0 and the temperature T .
The case of physical interest is the first one. Since in core of tokamaks, we have a large temperature T which implies a constant c T large, and a small colisionnality, in other words a small ν.
4 Existence of solutions.

:existence
In this section, we prove the existence theorem thm:existence 1.1. The proof will use the following notation and a preliminary lemma. A priori estimates of the Lemma where For each n > 0, we also introduce an approximation of the potential Φ n defined for any f ∈ L 2 m by

it holds that
for some constantc n . Substracting the equation satisfied by g 2 from the equation satisfied by g 1 and integrating over Ω leads to And so . Hence there is a unique fixed point of the map S on [0, T 1 ] for T 1 small enough. The bounds used for defining T 1 being independent of T 1 , a unique solution of the problem can be determined globally in time by iteration. The fact that this unique solution satisfy the a-priori estimates of the next section is clear since these estimates only depends on the bound satisfied by Φ and not on its precise form.

Proof of Theorem
The sequence (f n ) is compact in L 2 loc,u ((0, T ) × Ω). Indeed, it is bounded in L ∞ (0, T ; L 2 u (Ω U )) ∩ L 2 (0, T ; H 1 u (Ω U )) for any bounded subset U > 0 (Recall of Ω U = T 2 × R + ). It follows from the interpolation theory that (f n ) is bounded in L 10 3 u ((0, T ) × T 2 × U ). Together with the boundedness of (∇ x Φ n (f n )) in L 2 loc,u ((0, T ) × Ω), this implies that ( ∂fn ∂t ) is bounded in W   8], it holds that (f n ) is compact in L 2 u ((0, T ) × Ω U ), so converges up to a subsequence to some function f in L 2 u . It remains to pass to the limit when n → +∞ in the weak formulation satisfied by f n . A weak form of ( 0 uν∇ x f n · ∇ x α + ∂ u f n ∂ u α + βu 2 f n ∂ u α dxduds . (4.2) eq:weakform The passage to the limit in ( eq:weakform 4.2) when n → +∞ can be performed if This holds since (f n ) (resp. (∇ x (J 0 u (Φ n )(f n ))) strongly (resp. weakly) converges to f (resp. ∇ x (J 0 u (Φ(f ))) in L 2 loc,u ((0, T ) × Ω). And since the f n satisfy all the a priori bound, the limit f also satisfies them.

3.4) and the Jensen inequality
Integrating in u with the weight um, it holds using again Jensen inequality that 3.15) (we emphasize that the form used here is slighty different because we are in here in the 1D setting for the variable u) τ * 0 ∇f (t) 2 2,u 0 dt ≤ 2 A Appendix: Some controls on the Bessel function of zero-th order.

App:A
The first Bessel function J 0 is much used in this paper. Indeed, in Fourier space, J 0 1 that appears in the definition of the gyroaverage of the electric field, is the multiplication by J 0 . Some properties of the function J 0 are given in Watson [10]. In this appendix, some bounds on J 0 and its derivative are proven. The bound by ( √ 2k) − 1 2 is obtained as follows. J 0 is solution of the ordinary differential equation The new unknown u = √ kJ 0 is solution to There are no exact initial conditions for u. However, Moreover the series expansion of J 0 near k = 0, and its alternate character if k ≤ 2 imply that u 2 (k) ≥ k − k 3 2 (valid for k ≤ √ 2). Using, the inequality and the behavior of H near 0, we get if 0 < k 0 ≤ k ≤ √ 2 since the first line is satisfied for any k 0 > 0. Therefore, u 2 (k) ≤ k k 2 + 2 4k 2 + 1 , k ≤ √ 2 .