Diffusion limit of kinetic equations for multiple species charged particles
aa r X i v : . [ m a t h . A P ] J u l Diffusion limit of kinetic equations for multiple species chargedparticles
Hao Wu ∗ , Tai-Chia Lin † , and Chun Liu ‡ August 25, 2018
Abstract
In ionic solutions, there are multi-species charged particles (ions) with different proper-ties like mass, charge etc. Macroscopic continuum models like the Poisson–Nernst–Planck(PNP) systems have been extensively used to describe the transport and distribution of ionicspecies in the solvent. Starting from the kinetic theory for the ion transport, we study aVlasov–Poisson–Fokker–Planck (VPFP) system in a bounded domain with reflection bound-ary conditions for charge distributions and prove that the global renormalized solutions ofthe VPFP system converge to the global weak solutions of the PNP system, as the small pa-rameter related to the scaled thermal velocity and mean free path tends to zero. Our resultsmay justify the PNP system as a macroscopic model for the transport of multi-species ionsin dilute solutions.
Keywords : Ionic solutions, kinetic equation, diffusion limit, renormalized solution.
AMS Subject Classification : 35Q99, 35B25, 45K05, 35J05.
The transport of ions in different biological environments is very important in our life andit has attracted more and more attentions recently [25, 36, 37, 56]. In biological problems, theionic solutions usually consist of charged particles (ions) like sodium Na + , potassium K + , cal-cium Ca and chloride Cl − etc, which have different but comparable masses, charge valenciesand sizes. These differences have dramatic effects on the dynamics of multi-species ions whichproduce the functions of cells in biological system, e.g., the ion channels. To study the dynam-ics of multi-species ions, molecular dynamics simulations (MD) using microscopic models (fromNewton’s laws) to describe charge particle trajectories are popular and useful but expensive ∗ School of Mathematical Sciences and Shanghai Key Laboratory for Contemporary Applied Mathematics,Fudan University, 200433 Shanghai, China, Email: [email protected] † Department of Mathematics, National Taiwan University, No.1, Sec. 4, Roosevelt Road, Taipei 106, Taiwan,Email: [email protected] ‡ Department of Mathematics, Penn State University, State College, PA 16802, Email: [email protected] ∂ t c i = ∇ · J i ,J i = d i (cid:18) ∇ c i + q i k B T c i ∇ φ (cid:19) , −∇ · ( ǫ ∇ φ ) = P Ni =1 q i c i + D ( x ) . (1.1)Here, c i ( i = 1 , , ..., N ) stand for densities of charged particles in the ionic solution and φ is theself-consistent electric potential. Besides, q i are the (positive or negative) charges of particles, J i are the ionic flux densities, d i are their diffusion coefficients, ǫ is the dielectric coefficient, k B isthe Boltzmann constant, T is the temperature and D ( x ) is the permanent charge density in thedomain. The PNP system (1.1) is also one of the fundamental macroscopic models in the studyof transport of carriers in semiconductors, see, e.g., [31,43,47,48]. Concerning the mathematicalanalysis, the initial value problem and the initial boundary value problem of the PNP systemhave been extensively studied in the literature, we refer to [2,4,5,30–32,42,48] and the referencecited therein.The PNP system (1.1) provides a continuum description of the evolution of charged particlesvia macroscopic (averaged) quantities, e.g., the particle density, the current density etc., whichhave cheaper costs for numerics. Such continuum models can be (formally) derived from kineticmodels by coarse graining methods, like the moment method, the Hilbert expansion method andso on [42, 48, 50]. Although many results for the PNP systems have been obtained, it seems thatnone of them can reveal basic principles like gating and selectivity of ion channels. Recently, newPNP type systems have been derived and the selectivity of ion channels have been simulatedsuccessfully [27, 38–40, 46, 63]. In order to justify these continuum models, here we develop thekinetic theory for the PNP system like (1.1) as the first step work. Our goal in the presentpaper is to rigourously justify the PNP system for dilute ionic solutions consisting of multi-species charged particles, by studying the diffusion limit of a suitable kinetic system. We willcontinue to study the kinetic theory for those new PNP type systems as in [27, 38–40, 46, 63] inthe near future.In this paper, we consider the case that the motion of multi-species charged particles isgoverned by the electrostatic force coming from their (self-consistent) Coulomb interaction. Wealso assume that the momentum of charged particles with collision is small and ignorable. Thenthe collision term in the kinetic equation may be approximated by the Fokker–Planck operatorthat describes the Brownian force [18], and the resulting kinetic system becomes the Vlasov–Poisson–Fokker–Planck (VPFP) system as follows: ∂ t f i + v · ∇ x f i − z i qm i ∇ x φ · ∇ v f i = 1 τ i L iF P ( f i ) , (1.2)2 ǫ ∆ x φ = q N X i =1 z i Z R d f i dv + D ( x ) ! , (1.3)where L iF P ( i = 1 , ..., N ) are the Fokker–Planck operators such that L iF P ( f i ) = ∇ v · ( vf i + θ i ∇ v f i ) . Here, ǫ > q > i =1 , ..., N , the state of each species is given by a distribution function f i ( t, x, v ) ≥
0, i.e., aprobability density in the ( x, v )-phase space at time t ( f i dxdv is the number of the i -th speciescharged particles at time t located at a volume element dx about the position x and havingvelocities in a volume dv about the value v ). Besides, z i ∈ Z are the valencies for the N -speciescharged particles, m i are the masses, τ i are relaxation time due to collisions of the particleswith the thermal bath, √ θ i are the thermal velocities given by √ θ i = q k B T b m − i and T b is thetemperature of the thermal bath.In plasma physics, the VPFP system (1.2)–(1.3) with N = 1 (i.e., the single species case) isreasonable because the mass ratio between the ions and electrons is huge, only the evolution ofthe electrons is described in terms of a distribution function in the resulting system, and those‘heavy’ ions are supposed to be static. For such a case, the existence and uniqueness of solutionsto the initial value problem or the initial boundary value problem of the VPFP system havebeen investigated in the literature. We refer to [7, 58, 61] for results on the classical solutionsand to [8, 10, 11, 60] for weak solutions and their regularity. Concerning the long-time behaviorof the VPFP system, we refer to [6, 9, 12]. Instead of the single species case with N = 1, herewe study the system (1.2)–(1.3) with N ≥ L be the characteristic length. We denote by N the characteristic value for theconcentration of particles and by Φ the characteristic variation of the electric potential over L . Since we have to treat mutiple species of charged particles that have different masses andcharges, it is convenient to introduce a ‘reference particle’ with mass m ref , electric charge z ref q (with z ref = 1), relaxation time τ ref and thermal velocity θ ref . The microscopic variation aswell as the drift velocity for the reference particle are given by V ref = p θ ref , U ref = τ ref qm ref Φ L ,respectively. Choosing the following scaling (with respect to the reference particle) t → T t ′ , x → Lx ′ , v → V ref v ′ , T = LU ref and the change of unknowns f i ( t, x, v ) = N V − dref f ′ i ( t ′ , x ′ , v ′ ), φ ( t, x, v ) = Φ φ ′ ( t ′ , x ′ , v ′ ), D ( x ) = N D ′ ( x ′ ), we obtain the rescaled VPFP equations (drop theprime for simplicity): ∂ t f i + νv · ∇ x f i − κ i z i ε ∇ x φ · ∇ v f i = ζ i νε ∇ v · ( vf i + κ i ∇ v f i ) , i = 1 , ..., N, − ̟ ∆ x φ = N X i =1 z i Z R d f i dv + D ( x ) , ν (the ‘scaled’ thermal velocity), ε (the ‘scaled’ thermalmean free path), ̟ and the ratios κ i , ζ i are given by ν = V ref U ref , ε = τ ref V ref L , ̟ = ǫ Φ qN L , κ i = m ref m i , ζ i = τ ref τ i . The case we are interested in this paper is called the low field limit (or the parabolic limit),which means that the drift velocity is small comparing with the thermal velocity, while thethermal velocity is small comparing to the relaxation velocity, and the two ratios have the sameorder of magnitude (cf. [1, 29, 33, 57]): ν ≃ ε − and ε << . For ε >
0, taking ν = ε − (just for the sake of simplicity), we arrive at the rescaled VPFPsystem under low field scaling, which will be investigated in the remaining part of this paper: ∂ t f εi + 1 ε v · ∇ x f εi − κ i z i ε ∇ x φ ε · ∇ v f εi = ζ i ε L iF P ( f εi ) , (1.4) − ̟ ∆ x φ ε = N X i =1 z i Z R d f εi ( t, x, v ) dv + D ( x ) , (1.5)where the rescaled Fokker–Planck operators are given by L iF P ( f εi ) = ∇ v · ( vf εi + κ i ∇ v f εi ) . (1.6)We recall that the diffusion limit of the VPFP system has been studied extensively in theliterature (cf. [29, 33, 34, 55, 57] and the references therein). In [29, 33, 57], the authors studiedthe low field limit and proved the convergence of suitable solutions to the single species VPFPsystem towards a solution to the drift–diffusion–Poisson model in the whole space. In [57], undera suitable regularity assumption on the initial data, the convergence result was obtained globallyin time in two dimensions and locally in time for the three dimensional case. Later, the authorproved in [33] a global convergence result in the two dimensional case, without any restrictionon the time interval and the assumptions on the initial data were weakened with bounds onlyon the associated entropy and energy. Quite recently, in [29] the authors established a globalconvergence result, without any restriction on the time interval or on the spatial dimensions, byworking with the renormalized solutions (or free energy solutions, cf. [22, 24]). As pointed outin [29], the notion of renormalized solutions is natural for the problem, because the free energyof the VPFP system seems to be the only quantity that is uniformly bounded with respect tothe small parameter ε (i.e., the ‘scaled’ mean free path). Even one works with more regularinitial data such that the solutions can be defined in the usual weak sense without the needof renormalizing, one still has to use renormalization techniques to pass to the limit as ε → ∇ x φ · ∇ v f , where the maindifficulty comes from (we refer to [29] for more details).4n this paper, we rigorously prove that for the multi-species case, the VPFP system (1.4)–(1.5) converges to a rescaled PNP system as ε tends to zero in the low field limit. We generalizedthe techniques introduced in the previous works [29,49,53], to the case involving multiple speciesof charged particles in a bounded region with reflection boundary conditions [10, 15, 53]. Thespecific boundary conditions recover the classical no-flux boundary conditions of the PNP sys-tem. Different from the single species case in the literature, the previous arguments have to bemodified in order to deal with the nonlocal interactions between different species of particlesthrough the Poisson equation for the electric potential φ . Besides, in order to deal with theintegrals on the boundary, we shall make use of the Darroz`es–Guiraud information [21], whichhelps to obtain the energy dissipation. Finally, effects of different but comparable quantities likemasses and valencies of the charged particles will become obvious in our mathematical analysis.Our results support the PNP system (1.1) as a suitable model for multi-species chargedparticles in dilute solution. As we mentioned before, several variants of the PNP system (1.1)have recently been derived by using the energetic variational approaches [41] to model importantphysical ingredients such as size (steric) effects for non-diluted solutions (cf. e.g., [27, 38–40, 46,63]) that are crucial in the study of the selectivity of ion channels in cell membranes [19, 36,45]. The total energy for these modified PNP systems consists of the entropic energy inducedby the Brownian motion of ions, the electrostatic potential energy representing the coulombinteraction between the charged ions, and in particular, the repulsive potential energy causedby the excluded volume effect (e.g., the singular Lennard–Jones potential). Our result can beviewed as a starting point for the further investigation on the case of crowded ions. It would beinteresting to study the diffusion-limit of suitable kinetic systems to obtain the modified PNPsystems [27, 39].The remaining part of this paper is organized as follows. In Section 2, we present thedefinition of renormalized solutions and state the main result on the diffusion limit of the VFFPsystem (1.4)–(1.5) (Theorem 2.1). In Section 3, after deriving the energy dissipation of the VFFPsystem in bounded domain (Proposition 3.1), which yields the necessary uniform estimates(Lemmas 3.1, 3.2, 3.3), we proceed to prove our main result by using the renormalizationtechniques. Let Ω ⊂ R d ( d ≥
2) be a sufficiently smooth bounded domain. For instance, the outwardunit normal vector n ( x ) at x ∈ ∂ Ω satisfies n ∈ W , ∞ (Ω , R d ). The Lebesgue surface measureon ∂ Ω will be denoted by dS .Then we introduce the boundary conditions for the distribution functions. As in Cercignani’swork [14–16] (see also [6, 53]), we define the sets of outgoing (Σ x + ) and incoming (Σ x − ) velocitiesat point x ∈ ∂ Ω such that Σ x ± := { v ∈ R d : ± v · n ( x ) > } and denote the boundary sets5 ± = { ( x, v ) : x ∈ Ω , v ∈ Σ x ± } . Let γh be the trace of function h and γ ± h = (0 , + ∞ ) × Σ ± γh .Reflection boundary conditions for the kinetic equations take the form of integral (balance)relations between the densities of the particles on the outgoing and incoming velocity subsetsof the boundary ∂ Ω at a given time [14–16]. For instance, given x ∈ ∂ Ω and t >
0, we have(cf. [6]): γ − f ( t, x, v ) = Z Σ x + R ( t, x ; v, v ∗ ) γ + f ( t, x, v ∗ ) dv ∗ , v ∈ Σ x − , (2.1)where R represents the probability that a particle with velocity v ∗ at time t striking the boundaryon x reemerges at the same instant and location with velocity v . If we consider v ′ = − v forany v ∈ Σ x − and take R ( t, x ; v, v ∗ ) = δ v ′ being the Dirac measure centered at v ∗ = v ′ , thenwe have γ − f ( t, x, v ) = γ − f ( t, x, − v ) on Σ − , which is the classical (local) inverse reflectionboundary condition. Similarly, if we take v ′ = v − v · n ( x )) n ( x ), then we arrive at theclassical (local) specular reflection boundary condition, see [6, 10]. We refer to [6] for possibleminimal assumptions on R such that (2.1) is well-defined, i.e., R is nonnegative and it verifiesthe normalization condition as well as the reciprocity principle. Detailed discussions on theboundary conditions can be found in [14–16].Here, we are more interested in the so-called diffuse reflection according to a Maxwellianwith temperature of the thermal bath, which is nonlocal. Denote by M i ( v ) the Maxwellians forcharged particles M i ( v ) = 1(2 π ) d − κ d +12 i e − κi | v | , i = 1 , ..., N. (2.2)We note that M i are the zeros of the rescaled Fokker–Planck operators L iF P given in (1.6), i.e., L iF P ( M i ) = 0, ( i = 1 , ..., N ). Then we can choose a special form of R in (2.1) and propose thefollowing boundary conditions for the distribution functions (cf. [15,53]), which are special casesof the so-called Maxwell boundary condition [50, 53]: for given x ∈ ∂ Ω and t > γ − f εi = M i ( v ) R v · n ( x ) < | v · n ( x ) | M i ( v ) dv Z v ∗ · n ( x ) > ( γ + f εi ) v ∗ · n ( x ) dv ∗ , on Σ x − . (2.3)Besides, for the electric potential φ ε , we simply impose the zero-outward electric field conditionsuch that ∇ x φ ε · n = 0 , on ∂ Ω . (2.4)In summary, below we will consider the rescaled VPFP system (1.4)–(1.5) on (0 , T ) × Ω × R d subject to boundary conditions (2.3)–(2.4) and the initial data (depending on the parameter ε ): f εi ( t, x, v ) | t =0 = f εi ( x, v ) . (2.5)We remark that the boundary conditions (2.3) allow us to preserve mass conservation and obtainproper energy and entropy balance laws of the VPFP system (1.4)–(1.5). Denote by n εi ( t, x ) = Z R d f εi ( t, x, v ) dv and J εi = 1 ε Z R d vf εi dv, (2.6)6he densities as well as the current densities associated to the distribution functions for the i -thspecies, respectively. Multiplying (2.3) by v · n ( x ) and integrating over Σ x − , we easily deduce the(macroscopic) boundary conditions for the fluxes such that J εi · n = 0 , on ∂ Ω , (2.7)which imply that all the particles that reach the boundary are reflected (no particle goes outnor enters in the domain Ω) and thus the mass R Ω n εi ( t, x ) dx is conserved for all time. On theother hand, in order to uniquely determine the solution φ ε to the Poisson equation (1.5) withhomogeneous Neumann boundary condition (2.4), we require the global neutrality condition N X i =1 z i Z Ω Z R d f εi dvdx + Z Ω D ( x ) dx = 0 . (2.8)and the zero-mean constraint R Ω φ ε dx = 0. We first introduce the definition of renormalized solutions in the spirit of [29, 49]:
Definition 2.1.
The set ( f εi , φ ε ) ∈ L ∞ (0 , T ; ( L (Ω × R d )) N × H (Ω)) is a renormalized solutionto the VPFP system (1.4) – (1.5) with initial and boundary conditions (2.3) – (2.5) , if(1) For all functions β i ∈ C ( R ) , i = 1 , ..., N satisfying | β i ( s ) | ≤ C ( s + 1) , | β ′ i ( s ) | ≤ C (1 + s ) − , | β ′′ i ( s ) | ≤ C (1 + s ) − , ∀ s ≥ , the set ( β i ( f εi ) , φ ε ) is a weak solution to the system ε∂ t β i ( f εi ) + v · ∇ x β i ( f εi ) − κ i z i ∇ x φ ε · ∇ v β i ( f εi ) = ζ i ε L iF P ( f εi ) β ′ i ( f εi ) , (2.9) − ̟ ∆ x φ ε = N X i =1 z i Z R d f εi dv + D ( x ) , (2.10) with initial data β i ( f ε ) | t =0 = β i ( f εi ) (2.11) and boundary conditions γ − β i ( f εi ) = M i ( v ) R v · n ( x ) < | v · n ( x ) | M i ( v ) dv Z v ∗ · n ( x ) > γ + β i ( f εi ) v ∗ · n ( x ) dv ∗ , (2.12) ∇ x φ ε · n = 0 . (2.13) (2) For any λ > , the functions θ iε,λ = ( f εi + λ f M i ) satisfy ε∂ t θ iε,λ + v · ∇ x θ iε,λ − κ i z i ∇ v · ( ∇ x φ ε θ iε,λ ) = ζ i εθ iε,λ L iF P ( f εi ) + z i λ f M i θ iε,λ v · ∇ x φ ε , (2.14)7 here f M i are the normalized Maxwellians (comparing with (2.2) ) f M i ( v ) = (cid:16) κ i π (cid:17) M i ( v ) such that Z R d f M i ( v ) dv = 1 , i = 1 , ..., N. (2.15) Remark 2.1.
Due to the regularity of renormalized functions β i , the corresponding boundaryconditions (2.12) for the renormalized distribution functions make sense. We refer to [3, 17, 59](see also [6, 53]) for more detailed discussions about the traces of distribution functions on theboundary. Next, we consider the rescaled version of the PNP system (1.1): ∂ t n i + ∇ x · J i = 0 , (2.16) − ̟ ∆ x φ = N X i =1 z i n i + D ( x ) , (2.17)with density currents given by J i = − ζ i ∇ x n i − z i ζ i n i ∇ x φ (2.18)and subject to the following boundary conditions and initial conditions: J i · n = ∇ x φ · n = 0 , on (0 , T ) × ∂ Ω , (2.19) n i | t =0 = n i , in Ω . (2.20)Moreover, we require that Z Ω φdx = 0 and Z Ω N X i =1 z i n i + D ( x ) ! dx = 0 . Then we introduce the weak formulation of the PNP system (2.16)–(2.20).
Definition 2.2.
We say that the set ( n i , φ ) is a weak solution to the initial boundary valueproblem of the PNP system (2.16) – (2.20) , if n i ∈ L ∞ (0 , T ; L log L (Ω)) , √ n i ∈ L (0 , T ; H (Ω)) ,∂ t n i ∈ L (0 , T ; W − , (Ω)) , φ ∈ L (0 , T ; H (Ω)) , where the function space L log L (Ω) is given by L log L (Ω) := (cid:26) n : n ≥ , Z Ω n (1 + | log n | ) dx < + ∞ (cid:27) and the PNP system (2.16) – (2.17) is satisfied in the weak sense: for any u ∈ C ∞ ([0 , T ]; C ∞ (Ω)) , ψ ∈ L (0 , T ; ( H (Ω)) ′ ) , Z Ω n i ( t, · ) u ( t, · ) dx − Z Ω n i u (0 , · ) dx Z t Z Ω n i ∂ t udxdτ − ζ i Z t Z Ω ( ∇ x n i + z i n i ∇ x φ ) · ∇ x udxdτ, t ∈ [0 , T ] ,̟ Z T Z Ω ∇ x φ · ∇ x ψdxdt = Z T Z Ω N X i =1 z i n i + D ( x ) ! ψdxdt. Moreover, the weak solution ( n i , φ ) satisfies the following energy inequality e ( t ) + N X i =1 Z t Z Ω ζ i n i (cid:12)(cid:12)(cid:12) ∇ (cid:16) ln n i + z i φ (cid:17)(cid:12)(cid:12)(cid:12) dxdt ≤ e (0) , t ∈ [0 , T ] , with e ( t ) := Z Ω N X i =1 n i ln n i + ̟ |∇ φ | ! dx. Now we are in a position to state the main result of this paper.
Theorem 2.1.
Let the background charge D be independent of time and satisfy D ( x ) ∈ L ∞ (Ω) .We assume that the initial data ( f εi , φ ε ) satisfy the following assumptions f εi ≥ , Z Ω Z R d f εi (1 + | v | + | log f εi | ) dvdx ≤ C , (2.21) k φ ε k H (Ω) ≤ C , Z Ω φ ε dx = 0 , (2.22) for some constant C > independent of the parameter ε , and the global neutrality conditionholds N X i =1 z i Z Ω Z R d f εi dvdx + Z Ω D ( x ) dx = 0 , ∀ ε > . (2.23) Let ( f εi , φ ε ) be a free energy (renormalized) solution of the VPFP system (1.4) – (1.5) with cor-responding initial and boundary conditions (2.3) – (2.5) (cf. Definition 2.1). Then, as ε tends tozero, up to a subsequence if necessary, we have the strong convergence results f εi ( t, x, v ) → n i ( t, x ) f M i ( v ) in L (0 , T ; L (Ω × R d )) , (2.24) φ ε ( t, x ) → φ ( t, x ) in L (0 , T ; W ,p (Ω)) , ≤ p < . (2.25) Moreover, n εi strongly converge in L (0 , T ; L (Ω)) towards n i and ( n i , φ ) is a weak solution tothe PNP system (2.16) – (2.20) (cf. Definition 2.2) with initial data n i | t =0 = n i = R R d f i dv ,such that f i are the weak limits of f εi . Remark 2.2.
We would like to mention that the PNP system (2.16) – (2.20) can also be derivedfrom diffusion limits of other types of kinetic equations, e.g., the Boltzmann–Poisson system.We refer to [49] for the one species case and we believe that their argument can also be extendedto the multi-species case. Remark 2.3.
We remark that different types of scalings can be chosen for the VPFP system.For instance, if we assume that the drift and thermal velocities are comparable, but both are mall comparing with the relaxation velocity, e.g., ν = O (1) and ε << , then we arrive at adifferent rescaled VPFP system ∂ t f εi + v · ∇ x f εi − κ i z i ε ∇ x φ ε · ∇ v f εi = ζ i ε L iF P ( f εi ) , − ̟ ∆ x φ ε = n X i =1 z i Z R d f εi dv + D ( x ) . This is usually called drift-collision balance scaling or high field scaling in the literature. Takingthe hydrodynamic limit as ε → (the high field limit or the hyperbolic limit), the above VPFPsystem will lead to a first-order hyperbolic system for the density of particles coupled with thePoisson equation, cf. e.g., [1, 11, 34, 55]. The free energy of the VPFP system (1.4)–(1.5) is defined as follows E ( t ) = N X i =1 Z Ω Z R d (cid:18) | v | κ i f εi + H ( f εi ) (cid:19) dvdx + ̟ Z Ω |∇ x φ ε | dx, (3.1)where the function H takes the form H ( s ) = s log s for s ≥
0. The entropy productions of theVPFP system are given by D i ( w ) = Z Ω Z R d ( v √ w + 2 κ i ∇ v √ w ) dvdx = 4 Z Ω Z R d (cid:12)(cid:12)(cid:12)(cid:12)(cid:12) ∇ v q we κi | v | (cid:12)(cid:12)(cid:12)(cid:12)(cid:12) e − κi | v | dvdx, i = 1 , ..., N. (3.2)Moreover, we introduce the Darroz`es–Guiraud information on the boundary (cf. e.g., [21]) suchthat I i ( w ) = Z Σ x + H ( w ) dµ ix − H Z Σ x + wdµ ix ! , i = 1 , ..., N, where dµ ix ( v ) = M i ( v ) | v · n ( x ) | dv are probability measures on Σ x ± by the particular choice of thenormalized Maxwellians M i (cf. (2.2)).First, we derive the energy dissipation property of the VPFP system (1.4)–(1.5) with initialand boundary conditions (2.3)–(2.5). Proposition 3.1 (Energy dissipation) . The renormalized solution of the VPFP system (1.4) – (1.5) with described initial data and boundary conditions satisfies ∂ t n εi + ∇ x · J εi = 0 , (3.3)10 here n εi and J εi are given in (2.6) . Moreover, the following dissipative energy inequality holds E ( t ) + 1 ε N X i =1 ζ i κ i Z t D i ( f εi ) ds + 1 ε N X i =1 Z t Z ∂ Ω I i (cid:18) γ + f εi M i ( v ) (cid:19) dSds ≤ E (0) , ∀ t ≥ . (3.4) Proof.
We just present a formal calculation which leads to (3.4). For i = 1 , ..., N , multiplyingthe i -th equation in (1.4) of the VPFP system by | v | and integrating the result with respectto v and x , we get ddt Z Ω Z R d | v | f εi dvdx + Z Ω Z R d ε | v | v · ∇ x f εi dvdx − Z Ω Z R d κ i z i ε | v | ∇ x φ ε · ∇ v f εi dvdx = Z Ω Z R d ζ i ε | v | L iF P ( f εi ) dvdx, integrating by parts, we see that Z Ω Z R d ε | v | v · ∇ x f εi dvdx = 12 ε Z ∂ Ω Z R d ( v · n ) | v | γf εi dvdS, − Z Ω Z R d ε | v | ∇ x φ ε · ∇ v f εi dvdx = 1 ε Z Ω Z R d ( v · ∇ x φ ε ) f εi dvdx = − Z Ω φ ε ∇ x · J εi dx + Z ∂ Ω γφ ε J εi · n dS = Z Ω φ ε ∂ t n εi dx, and Z Ω Z R d ε | v | L iF P ( f εi ) dvdx = − ε Z Ω Z R d ( vf εi + κ i ∇ v f εi ) · vdvdx. As a result, for i = 1 , ..., N we obtain that ddt Z Ω Z R d κ i | v | f εi dvdx + z i Z Ω φ ε ∂ t n εi dx = − κ i ε Z ∂ Ω Z R d ( v · n ) | v | γf εi dvdS − ζ i κ i ε Z Ω Z R d ( vf εi + κ i ∇ v f εi ) · vdvdx. (3.5)Next, multiplying the i -th equation (1.4) of the VPFP system by log f εi and integrating theresult with respect to v and x , we get ddt Z Ω Z R d H ( f εi ) dvdx + Z Ω Z R d ε ( v · ∇ x f εi ) log f εi dvdx − Z Ω Z R d κ i z i ε ( ∇ x φ ε · ∇ v f εi ) log f εi dvdx Z Ω Z R d ζ i ε L iF P ( f εi ) log f εi dvdx, after integrating by parts, we see that Z Ω Z R d ε ( v · ∇ x f εi ) log f εi dvdx = − ε Z Ω Z R d ( v · ∇ x f εi ) dvdx + 1 ε Z ∂ Ω Z R d ( v · n ) γf εi log γf εi dvdS = − Z ∂ Ω J εi · n dS + 1 ε Z ∂ Ω Z R d ( v · n ) γf εi log γf εi dvdS = 1 ε Z ∂ Ω Z R d ( v · n ) γf εi log γf εi dvdS, − Z Ω Z R d ε ( ∇ x φ ε · ∇ v f εi ) log f εi dvdx = 0 , and Z Ω Z R d ζ i ε L iF P ( f εi ) log f εi dvdx = − ζ i ε Z Ω Z R d ( vf εi + κ i ∇ v f εi ) · ∇ v f εi f εi dvdx. As a consequence, we get ddt Z Ω Z R d H ( f εi ) dvdx = − ε Z ∂ Ω Z R d ( v · n ) γf εi log γf εi dvdS − ζ i ε Z Ω Z R d ( vf εi + κ i ∇ v f εi ) · ∇ v f εi f εi dvdx. (3.6)Moreover, we see that for f εi Z Ω Z R d ( vf εi + κ i ∇ v f εi ) · (cid:18) v + κ i ∇ v f εi f εi (cid:19) dvdx = D i ( f εi ) . (3.7)Due to the Poisson equation (1.5), we have Z Ω φ ε ∂ t N X i =1 z i n εi ! dx = − ̟ Z Ω φ ε ∂ t ∆ x φ ε dx = ̟ ddt Z Ω |∇ x φ ε | dx. (3.8)Then we conclude from (3.5)–(3.8) that ddt E ( t ) + N X i =1 ζ i κ i ε D i ( f εi ) + 1 ε N X i =1 Z ∂ Ω Z R d ( v · n ) (cid:18) κ i | v | + log γf εi (cid:19) γf εi dvdS = 0 . Recall that dµ ix ( v ) = M i ( v ) | v · n ( x ) | dv are probability measures on Σ x ± (see the definitionof M i ( v ) (2.2) and (2.15)). Then for the boundary terms, we can apply the Darroz`es–Guiraudinequality [21], namely, thanks to (2.3), the convexity of H ( s ) = s log s and the Jensen inequalitywe deduce that (see also [53]) Z ∂ Ω Z R d ( v · n ) (cid:18) κ i | v | + log γf εi (cid:19) γf εi dvdS Z ∂ Ω Z Σ x + H (cid:18) γ + f εi M i ( v ) (cid:19) dµ ix dS − Z ∂ Ω Z Σ x − H (cid:18) γ − f εi M i ( v ) (cid:19) dµ ix dS = Z ∂ Ω Z Σ x + H (cid:18) γ + f εi M i ( v ) (cid:19) dµ ix dS − Z ∂ Ω H Z Σ x + γ + f εi M i ( v ) dµ ix ! dS = Z ∂ Ω I i (cid:18) γ + f εi M i ( v ) (cid:19) dS ≥ . As a consequence, ddt E ( t ) + 1 ε N X i =1 ζ i κ i D i ( f εi ) + 1 ε N X i =1 Z ∂ Ω I i (cid:18) γ + f εi M i ( v ) (cid:19) dS ≤ . (3.9)Integrating (3.9) with respect to time, we arrive at our conclusion (3.4).The energy dissipation (3.4) yields the following global estimates that are uniform in theparameter ε , which enable us to take the diffusion limit as ε → Lemma 3.1.
For any
T > , there exists a constant C depending on C , ζ i , κ i , ̟ , but inde-pendent of ε and t ∈ [0 , T ] such that Z Ω Z R d (1 + | v | + | log( f εi ) | ) f εi dvdx ≤ C, Z Ω |∇ x φ ε | dx ≤ C, ε Z t D i ( f εi ) ds ≤ C, ε Z t Z ∂ Ω I i (cid:18) γ + f εi M i ( v ) (cid:19) dSds ≤ C. The functions f εi are weakly relatively compact in L ((0 , T ) × Ω × R d ) and fulfill k∇ v p f εi k L ((0 ,T ) × Ω × R d ) ≤ C. Concerning the fluxes, we have k J εi ( t, · ) k L (Ω) ≤ ε D i ( f εi ) + 12 k f εi k L (Ω × R d ) . Proof.
The proof is similar to [29, Propositions 5.1, 5.2, 5.3], based on the energy inequality(3.4). Since we are now dealing with the bounded domain, we do not need to estimate termslike R Ω R R d | x | f εi dvdx as in [29]. The L weak compactness of f εi follows from the well-knownDunford–Pettis theorem.We recall that the initial boundary value problem of a full Vlasov–Poisson–Fokker–Planck–Boltzmann system (subject to more general reflection boundary conditions for the distribution13unction but only for one species of charged particles) has been studied in the recent paper [53].The author proved the existence of DiPerna–Lions renormalized solutions by using the approx-imation procedure in [51] with some crucial trace theorems previously introduced by the sameauthor for the Vlasov equations [52] and some new results concerning weak-weak convergence(the renormalized convergence and the biting L -weak convergence). For the current case withmultiple species of charged particles, the coupling between different species is somewhat weak,i.e., only via the Poisson equation. As a result, based on the energy dissipation property Propo-sition 3.1 and Lemma 3.1, we are able to prove the following existence result on renormalizedsolutions to the VPFP system (1.4)–(1.5), by adapting the argument in [53] (see also [9, 49, 51])with minor modifications. The details are thus omitted. Theorem 3.1 (Existence of renormalized solution) . Suppose that the assumptions (2.21) – (2.23) on the initial data are satisfied. For arbitrary but fixed ε > , the initial boundary value problemof the VPFP system (1.4) – (1.5) admits at least one (renormalized) solution ( f εi , φ ε ) in the senseof Definition 2.1, which satisfies Proposition 3.1. ε → The proof of Theorem 2.1 mainly follows the arguments in [29] for the VPFP system thatconcerns only one single species of particles in the whole space. However, for the present probleminvolving multiple species of charged particles, we need to modify the previous argument to dealwith nonlocal interactions between particles as well as the boundary conditions. In what follows,we state the essential steps and point out the possible differences in the proof.
Step 1. Strong convergence of the electric potential φ ε .Based on the uniform estimates in Lemma 3.1, it is straightforward to argue as [49, Propi-sition 3.3] to conclude that Lemma 3.2.
The renormalized solution ( f εi , φ ε ) satisfies the following properties:(1) for i = 1 , ..., N , n εi ( t, x ) = R R d f εi ( t, x, v ) dv are weakly relatively compact in L ((0 , T ) × Ω) ,(2) φ ε ( t, x ) is relatively compact in L (0 , T ; W ,p (Ω)) with ≤ p < . Therefore, the strong convergence of φ ε (2.25) (up to a subsequence) is a direct consequenceof Lemma 3.2. Step 2. Strong convergence of the charge densities n εi .Lemma 3.2 also implies the weak compactness of densities n εi . Indeed, we can show theconvergence of density functions in the strong sense. By using the definition of renormalizedsolutions (cf. Definition 2.1) and a velocity averaging lemma (cf. [49, Lemma 4.2], also [23]),we are able to obtain the compactness of the densities (cf. [29, Proposition 6.1]) such that thedensities n εi are relatively compact in L ((0 , T ) × Ω), namely, there exist n i ∈ L ((0 , T ) × Ω)and up to a subsequence if necessary, n εi → n i , in L ((0 , T ) × Ω) and a.e. as ε → . (3.10)14he above result and the simple inequality ( √ a − √ b ) ≤ | a − b | imply that p n εi → √ n i , in L ((0 , T ) × Ω) and a.e. as ε → . (3.11) Step 3. Strong convergence of the distribution functions f εi .We recall the classical Csiszar–Kullback inequality (cf. [20, Theorem 3.1, Section 4, pp. 314],see also [44]) that for all non-negative u ∈ L ( R d , dµ ) (where dµ is a probability measure) with R R d udµ = 1, it holds k u − k L ( R d ,dµ ) ≤ (cid:18)Z R d ( u log u − u + 1) dµ (cid:19) . Choose in the above inequality u = f εi n εi f M i ( v ) , dµ = f M i ( v ) dv, which easily implies (cid:18)Z Ω Z R d | f εi − n εi f M i ( v ) | dvdx (cid:19) ≤ (cid:18)Z Ω n εi dx (cid:19) Z Ω Z R d f εi log f εi n εi f M i ( v ) ! dvdx. (3.12)Next, we proceed to estimate the second factor in the righthand side of (3.12). Recalling thelogarithmic Sobolev inequality (cf. e.g., [35, Corollary 4.2]) Z R d | h ( v ′ ) | log | h ( v ′ ) | dµ ( v ′ ) ≤ Z R d |∇ v ′ h ( v ′ ) | dµ ( v ′ ) + k h ( v ′ ) k L ( R d ,dµ ( v ′ )) log k h ( v ′ ) k L ( R d ,dµ ( v ′ )) , where dµ ( v ′ ) is the Gauss measure dµ ( v ′ ) = (2 π ) − d e − | v ′| dv . Making the simple change ofvariable v ′ → v √ κ and denoting h κ ( v ) = h ( v ′ ), we have dµ ( v ′ ) = (cid:16) κ π (cid:17) d e − κ | v | dv := dµ κ ( v ) , k h ( v ′ ) k L ( R d ,dµ ( v ′ )) = k h κ ( v ) k L ( R d ,dµ κ ( v )) , which yields that Z R d | h κ ( v ) | log | h κ ( v ) | dµ κ ( v ) ≤ κ Z R d |∇ v h κ ( v ) | dµ κ ( v ) + k h κ ( v ) k L ( R d ,dµ κ ( v )) log k h κ ( v ) k L ( R d ,dµ κ ( v )) . In the above inequality, we set κ = κ i , h κ ( v ) = s f εi f M i ( v ) , dµ κ ( v ) = f M i ( v ) dv. k h κ k L ( R d ,dµ κ ( v )) = n εi , which yields Z Ω Z R d f εi log f εi n εi f M i ( v ) ! dvdx = Z Ω (cid:18) Z R d | h κ ( v ) | log | h κ ( v ) | dµ κ ( v ) − k h κ ( v ) k L ( R d ,dµ κ ( v )) log k h κ ( v ) k L ( R d ,dµ κ ( v )) (cid:19) dx ≤ κ i Z Ω Z R d (cid:12)(cid:12)(cid:12)(cid:12)(cid:12) ∇ v s f εi f M i ( v ) (cid:12)(cid:12)(cid:12)(cid:12)(cid:12) f M i ( v ) dvdx = κ i D i ( f εi ) . (3.13)As a consequence, we infer from the entropy dissipation in (3.4), the uniform estimates inLemma 3.1 and the estimates (3.12) and (3.13) that when ε → f εi − n εi f M i → , in L ((0 , T ) × Ω × R d ) and a.e.Combing the above results with the convergence result of n εi (3.10), we conclude that as ε → f εi → n i f M i , in L ((0 , T ) × Ω × R d ) and a.e. (3.14)Here and below, the convergence results are always understood to be up to a subsequence. Step 4. Weak convergence of the fluxes J εi .We introduce the auxiliary functions r εi = p f εi − q n εi f M i ( v ) ε q f M i ( v ) , i = 1 , ..., N. (3.15)In analogy to [49, Proposition 3.4] and [29, Proposition 5.5], we have Lemma 3.3.
For arbitrary
T > , the following uniform estimates hold Z T Z Ω Z R d (cid:16) | r εi | f M i + ε | r εi | | v | f M i + √ ε | r εi | | v | f M i (cid:17) dvdxdt ≤ C, where C is a constant that may depend on C , ζ i , κ i , ̟ , but independent of ε and t ∈ [0 , T ] . Using the expressions of r εi (cf. (3.15)), we have f εi = n εi f M i + 2 ε f M i p n εi r εi + ε | r εi | f M i . (3.16)Due to the simple facts R R d v f M i ( v ) dv = 0, it follows from (2.6), (3.11) and Lemma 3.3 that as ε → J εi = 2 p n εi Z R d r εi v f M i dv + √ ε Z R d √ ε | r εi | v f M i dv √ n i Z R d r i v f M i dv, weakly in L ((0 , T ) × Ω) , where r i are the weak limits of r εi , for i = 1 , ..., N .It remains to identify the limit function of J εi , which can be done by using a similar argumentas in [29, Proposition 7.2]. The strong convergence of f εi (see (3.14)) implies that for any fixed λ > θ iε,λ → q ( n i + λ ) f M i , as ε → . On the other hand, it follows from (3.16) that for any λ >
0, when ε →
0, we have ζ i ε L iF P ( f εi ) = ζ i L iF P (cid:16) f M i p n εi r εi + ε | r εi | f M i (cid:17) → ζ i √ n i L iF P ( r i f M i ) . As a consequence, in the renormalized formula (2.14), first for any fixed λ > ε → λ →
0, we obtain that (cid:16) ∇ x √ n i + z i ∇ x φ √ n i (cid:17) · v f M i = ζ i L iF P ( r i f M i ) , (3.17)where φ is the limit of φ ε (recall (2.25)).On the other hand, it follows from [29, Proposition 3.1] that χ j = − v j f M i ( i = 1 , ..., N , j = 1 , ..., d ) is the unique solution to the equation L iF P χ j = v j f M i in R ( L iF P ) ∩ D ( L iF P ), where L f M i ( R d ) = L ( R d ; f M − i dv ) ,R ( L iF P ) = (cid:26) f ∈ L f M i ( R d ) : Z R d f ( v ) dv = 0 (cid:27) ,D ( L iF P ) = n f ∈ L f M i ( R d ) : ∇ v · (cid:16) e − κi | v | ∇ v ( e κi | v | f ) (cid:17) ∈ L f M i ( R d ) o . Since − L iF P is a self-adjoint operator on L f M i ( R d ), using (3.17), we have J i = 2 √ n i Z R d r i v f M i dv = 2 √ n i Z R d ( r i f M i ) L iF P ( − v f M i ) f M − i dv = 2 √ n i Z R d L iF P ( r i f M i )( − v f M i ) f M − i dv = 2 ζ i √ n i Z R d h ( ∇ x √ n i + z i ∇ x φ √ n i ) · v f M i i ( − v f M i ) f M − i dv = − ζ i √ n i (cid:18)Z R d v ⊗ v f M i dv (cid:19) (cid:16) ∇ x √ n i + z i ∇ x φ √ n i (cid:17) = − ζ i √ n i (cid:16) ∇ x √ n i + z i ∇ x φ √ n i (cid:17) . where we use the fact that R R d v ⊗ v f M i dv = I . Therefore, we can see that as ε → J εi → J i := − ζ i √ n i (cid:16) ∇ x √ n i + z i ∇ x φ √ n i (cid:17) . (3.18)17n the distribution sense. Step 5. Passage to the limit in the PDE system .In order to recover the PNP system (2.16)–(2.20), we state a regularity result for the densityfunctions n i in the spirit of [49, Lemma 7.1] Lemma 3.4.
Let Ω be a smooth bounded and open set in R d . Assume n i are positive functionsbelonging to L ∞ (0 , T ; L (Ω)) and φ ∈ L (0 , T ; H (Ω)) that satisfy ∇ x √ n i + z i ∇ x φ √ n i = G i ∈ L (0 , T ; L (Ω)) , i = 1 , ..., N, (3.19) − ̟ ∆ x φ = N X i =1 z i n i + D ( x ) . Then we have √ n i ∈ L (0 , T ; H (Ω)) , N X i =1 z i n i ∈ L (0 , T ; L (Ω)) , ∇ x φ √ n i ∈ L (0 , T ; L (Ω)) . Proof.
As in [49, Corollary 3.2], we take β δ ( s ) = δ − β ( δs ) where β ∈ C ∞ ( R ) satisfying β ( s ) = s for − ≤ s ≤
1, 0 ≤ β ′ ( s ) ≤ s ∈ R and β ( s ) = 2 for | s | ≥
3. Then we renormalize theequations (3.19) for √ n i such that ∇ x β δ ( √ n i ) + z i ∇ x φβ ′ δ ( √ n i ) √ n i = G i β ′ δ ( √ n i ) ∈ L (0 , T ; L (Ω)) . (3.20)For any δ >
0, due to our choice of β and the given regularity for ∇ x φ , we have k∇ x φβ ′ δ ( √ n i ) √ n i k L (0 ,T ; L (Ω)) ≤ δ k∇ x φ k L (0 ,T ; L (Ω)) , which together with (3.20) implies that ∇ x β δ ( √ n i ) ∈ L (0 , T ; L (Ω)). Then we can take L norm on both sides of the equations (3.20), summing up with respect to i = 1 , ..., N , we have N X i =1 k∇ x β δ ( √ n i ) k L (0 ,T ; L (Ω)) + N X i =1 z i k∇ x φβ ′ δ ( √ n i ) √ n i k L (0 ,T ; L (Ω)) + N X i =1 Z T Z Ω (cid:2) z i β ′ δ ( √ n i ) √ n i ∇ x β δ ( √ n i ) (cid:3) · ∇ x φdxdt ≤ N X i =1 k G i k L (0 ,T ; L (Ω)) , where the right-hand side is independent of δ . For the crossing term on the left hand side, usingintegration by parts, we have N X i =1 Z T Z Ω (cid:2) z i β ′ δ ( √ n i ) √ n i ∇ x β δ ( √ n i ) (cid:3) · ∇ x φdxdt N X i =1 Z T Z Ω ∇ x h z i ˜ β δ ( √ n i ) i · ∇ x φdxdt = 1 ̟ Z T Z Ω N X i =1 z i ˜ β δ ( √ n i ) · N X i =1 z i n i + D ( x ) ! dxdt where ˜ β satisfies˜ β ( s ) = Z s τ β ′ ( τ ) dτ, ˜ β δ ( s ) = δ − ˜ β ( δs ) , ˜ β δ ( s ) → s , as δ → . Let δ →
0, we infer from the above estimates that N X i =1 k∇ x √ n i k L (0 ,T ; L (Ω)) + N X i =1 z i k∇ x φ √ n i k L (0 ,T ; L (Ω)) + 12 ̟ Z T Z Ω N X i =1 z i n i ! dxdt ≤ N X i =1 k G i k L (0 ,T ; L (Ω)) + 12 ̟ (cid:12)(cid:12)(cid:12)(cid:12)(cid:12)Z T Z Ω D ( x ) N X i =1 z i n i ! dxdt (cid:12)(cid:12)(cid:12)(cid:12)(cid:12) ≤ N X i =1 k G i k L (0 ,T ; L (Ω)) + 14 ̟ Z T Z Ω | D ( x ) | + N X i =1 z i n i ! dxdt, which easily yields the required regularity estimate. The lemma is proved.Finally, using the above regularity lemma and the convergence result (3.18), we are able towrite the currents J i as in (2.18). Then we can pass to the limit as ε → n i , φ ) satisfy the rescaled PNP system (2.16)–(2.20).The proof of Theorem 2.1 is complete. Acknowledgement
H. Wu was partially supported by NSF of China 11371098 and “Zhuo Xue” program ofFudan University. C. Liu was partially supported by NSF grants DMS-1109107, DMS-1216938and DMS-1159937.
References [1] A. Arnold, J.-A. Carrillo, I. Gamba and C.-W. Shu, Low and high field scaling limits forthe Vlasov– and Wigner–Poisson–Fokker–Planck systems, Transport Theory Statist. Phys., (2&3) (2001), 121–153. 192] A. Arnold, P. A. Markowich and G. Toscani, On large time asymptotics for drift–diffusion–Poisson systems, Transp. Th. Statist. Phys., (2000), 571–581.[3] R. Beals and V. Protopopescu, Abstract time dependent transport equations, J. Math.Anal. Appl., (1987), 370–405.[4] P. Biler and W. Hebisch and T. Nadzieja, The Debye system: existence and large timebehavior of solutions, Nonlinear Anal., (9) (1994), 1189–1209.[5] P. Biler and J. Dolbeault, Long time behavior of solutions of Nernst–Planck and Debye-H¨uckel drift–diffusion systems, Ann. Henri Poincar´e, (2000), 461–472.[6] L. Bonilla, J.-A. Carrillo and J. Soler, Asymptotic behavior of an initial-boundary valueproblem for the Vlasov–Poisson–Fokker–Planck system, SIAM J. Appl. Math., (5) (1997),1343–1372.[7] F. Bouchut, Existence and uniqueness of a global smooth solution for the Vlasov–Poisson–Fokker–Planck system in three dimensions, J. Funct. Anal., (1) (1993), 239–258.[8] F. Bouchut, Smoothing effect for the non-linear Vlasov–Poisson–Fokker–Planck system, J.Differ. Equs., (2) (1995), 225-238.[9] F. Bouchut and J. Dolbeault, On long asymptotics of the Vlasov–Fokker–Planck equationand of the Vlasov–Poisson–Fokker–Planck system with coulombic and newtonian potentials,Differential Integral Equations, (1995), 487–514.[10] J.-A. Carrillo, Global weak solutions for the initial-boundary value problems to the Vlasov–Poisson–Fokker–Planck system, Math. Meth. Appl. Sci., (1998), 907–938.[11] J.-A. Carrillo and J. Soler, On the initial value problem for the Vlasov–Poisson–Fokker–Planck system with initial data in L p spaces, Math. Methods Appl. Sci., (10) (1995),825–839.[12] J.-A. Carrillo, J. Soler and J.-L. Vazquez, Asymptotic behaviour and self-similarity forthe three dimensional Vlasov–Poisson–Fokker–Planck system, J. Funct. Anal., (1996),99–132.[13] C. Cercignani, The Boltzmann Equation and its Applications , Springer, Berlin, 1988.[14] C. Cercignani, Scattering kernels for gas/surface interaction, in Proceedings of the workshopon hypersonic flows for reentry problems, 1, INRIA, Antibes (1990), 9–29.[15] C. Cercignani, I. Gamba and C. Levermore, A drift-collision balance for a Boltzmann–Poisson system in bounded domains, SIAM J. Appl. Math., (6) (2001), 1932–1958.[16] C. Cercignani, M. Lampis and A. Lentati, A new scattering kernel in kinetic theory of gases,Transp. Theory Stat. Phys., (1995), 1319–1336.[17] Cessenat M., Th´eor`emes de trace pour des espaces de fonctions de la neutronique, CRAS, (3) (1985), 89–92. 2018] S. Chandrasekhar, Brownian motion, dynamical friction and stellar dynamics, Rev. Mod.Phys., (1949), 383–388.[19] R. Coalson and M. Kurnikova, Poisson–Nernst–Planck theory approach to the calculationof current through biological ion channels, IEEE transactions on nanobioscience, (2005),81–93.[20] I. Csiszar, Information-type measures of difference of probability distributions and indirectobservations, Stud. Sci. Math. Hung., (1967), 299–318.[21] J.-S. Darroz`es and J.-P. Guiraud, G´en´eralisation formelle du th´eor`eme H en pr´esence deparois, C.R.A.S. (Paris), A262 (1966), 1368–1371.[22] R. DiPerna and P.-L. Lions, Solutions globales d´equations du type Vlasov–Poisson, C.R.Acad. Sci. Paris Ser. I Math., (12) (1988), 655–658.[23] R.-J. Diperna, P.-L. Lions and Y. Meyer, L p regularity of velocity averages, Ann. Inst. H.Poincar´e Anal. Non Lin´eaire, (1991), 271–287.[24] J. Dolbeault, Free energy and solutions of the Vlasov–Poisson–Fokker–Planck system: ex-ternal potential and confinement (large time behavior and steady states), J. Math. PuresAppl., (9) (2) (1999), 121–157.[25] B. Eisenberg, Ionic channels in biological membranes: natural nanotubes, Acc. Chem. Res., (1998), 117–123.[26] B. Eisenberg, Multiple scales in the simulation of ion channels and proteins, J. PhysicalChemistry, (2010), 20719–20733.[27] B. Eisenberg, Y. Hyon and C. Liu, Energy variational analysis of ions in water and channels:field theory for primitive models of complex ionic fluids, J. Chemical Physics, (2010),104104.[28] B. Eisenberg and W. Liu, Poisson–Nernst–Planck systems for ion channels with permanentcharges, SIAM J. Math. Anal., (6) (2007), 1932–1966.[29] N. El Ghani and N. Masmoudi, Diffusion limit of the Vlasov–Poisson–Fokker–Planck sys-tem, Commun. Math. Sci., (2) (2010), 463–479.[30] W. Fang and K. Ito, On the time-dependent drift-diffusion model for semiconductors, J.Differential Equations, (1995), 245–280.[31] H. Gajewski, On existence, uniqueness, and asymptotic behavior of solutions of the basicequations for carrier transport in semiconductors, Z. Angew. Math. Mech., (1985), 101–108.[32] H. Gajewski and K. Gr¨oger, On the basic equations for carrier transport in semiconductors,J. Math. Anal. Appl., (1) (1986), 12–35.[33] T. Goudon, Hydrodynamic limit for the Vlasov–Poisson–Fokker–Planck system: analysisof the two-dimensional case, Math. Models Methods Appl. Sci., (5) (2005) 737–752.2134] T. Goudon, J. Nieto, F. Poupaud and J. Soler, Multidimensional high-field limit ofthe electro-static Vlasov–Poisson–Fokker–Planck system, J. Differetial Equations, (2)(2005), 418–442.[35] L. Gross, Logarithmic Sobolev inequalities, American J. Math., (4) (1975), 1061–1083.[36] B. Hille, Ion Channels of Excitable Membranes , 3rd ed. Sunderland, MA: Sinauer Asso-ciates, Inc., 2001.[37] A.-L. Hodgkin and A.-F. Huxley, A qualitative description of the membrane current andits application to conduction and excitation nerve, J. Physiol., (1952), 500–544.[38] T.-L. Horng, T.-C. Lin, C. Liu and B. Eisenberg, PNP equations with steric effects: amodel of ion flow through channels, J. Phys. Chem. B, (2012), 11422–11441.[39] Y. Hyon, B. Eisenberg and C. Liu, A Mathematical model for the hard sphere repulsion inionic solutions, Commun. Math. Sci., (2) (2011), 459–475.[40] Y. Hyon, J. Fonseca, B. Eisenberg and C. Liu, Energy variational approach to study chargeinversion (layering) near charged walls, Discrete Contin. Dyn. Syst. Ser. B, (8) (2012),2725–2743.[41] Y. Hyon, D.-Y. Kwak and C. Liu, Energetic variational approach in complex fluids: maxi-mum dissipation principle, Discrete Conti. Dynam. Sys., (4) (2010), 1291–1304.[42] A. J¨ungel, Quasi-hydrodynamic Semiconductor Equations , Progress in Nonlinear Differen-tial Equations and their Applications, , Birkhuser Verlag, Basel, 2001.[43] J.-W. Jerome, Analysis of Charge Transport – a Mathematical Study of SemiconductorDevices , Springer, 1996.[44] S. Kullback, A lower bound for discrimination information in terms of variation, IEEETrans. Infor. Theory, (1967), 126–127.[45] W. Kunz, Specific Ion Effects , World Scientific Publishing, Singapore, 2009.[46] T.-C. Lin and B. Eisenberg, A new approach to the Lennard–Jones potential and a newmodel: PNP-steric equations, Comm. Math. Sci., (1) (2014), 149–173.[47] P.-A. Markowich, The Stationary Semiconductor Device Equations , Computational Micro-electronics, Springer-Verlag, Vienna, 1986.[48] P.-A. Markowich, C.-A. Ringhofer and C. Schmeiser,
Semiconductor Equations , Springer,1990.[49] N. Masmoudi and M. Tayeb, Diffusion limit of a semiconductor Boltzmann–Poisson system,SIAM J. Math. Anal., (6) (2007), 1788–1807.[50] J. Maxwell, On stresses in rarefied gases arising from inequalities of temperature, Phil.Trans. Roy. Soc. London, (1879), Appendix 231–256.[51] S. Mischler, On the initial boundary value problem for the Vlasov–Poisson–Boltzmannsystem, Comm. Math. Phys., (2000), 447–466.2252] S. Mischler, On the trace problem for solutions of the Vlasov equation, Comm. PartialDifferential Equations, (2000), 1415–1443.[53] S. Mischler, Kinetic equations with Maxwell boundary conditions, Ann. Sci. ´Ec. Norm.Sup´er., (5) (2010), 719–760.[54] W. Nernst, Die elektromotorische Wirksamkeit der Ionen, Z. Phys. Chem., (1889), 129–181.[55] J. Nieto, F. Poupaud and J. Soler, High-field limit for the Vlasov–Poisson–Fokker–Plancksystem, Arch. Ration. Mech. Anal., (1) (2001), 29–59.[56] W. Nonner and B. Eisenberg, Ion permeation and glutamate residues linked by Poisson–Nernst–Planck theory in L-type calcium channels, Biophys. J., (1998), 1287–1305.[57] F. Poupaud and J. Soler, Parabolic limit and stability of the Vlasov–Fokker–Planck system,Math. Models Methods Appl. Sci., (7) (2000), 1027–1045.[58] G. Rein and J. Weckler, Generic global classical solutions of the Vlasov–Fokker–Planck–Poisson system in three dimensions, J. Differ. Equs., (1) (1992), 59–77.[59] S. Ukai, Solutions of the Boltzmann equation, in patterns and waves-qualitative analysis ofdifferential equations, Stud. Math. Appl., (1986), 37–96.[60] H.D. Victory, On the existence of global weak solutions for Vlasov–Poisson–Fokker–Plancksystems, J. Math. Anal. Appl., (2) (1991), 525–555.[61] H.D. Victory and B.P. O’Dwyer, On classical solutions of Vlasov–Poisson–Fokker–Plancksystems, Indiana Univ. Math. J., (1) (1990), 105–156.[62] L. Wan, S.-X. Xu, M. Liao, C. Liu and P. Sheng, Self-consistent approach to global chargeneutrality in electrokinetics: a surface potential trap model, Phys. Rev. X, (2014), 011042.[63] S.-X. Xu, P. Sheng and C. Liu, An energetic variational approach for ion transport, Com-mun. Math. Sci.,12