24P044-P
Sheath stability conditions for PIC simulation
PIC シミュレーションにおけるシース安定条件
Hiroto Matsuura,Naoki Inagaki
松浦 寛人, 稲垣直輝
Res. Org. for Univ.-Comm. Collabo., Osaka Pref. U.,
Gakuen-cho 1-2, Naka-ku, Sakai, Osaka 599-8570, Japan
大阪府立大学地域連携研究機構、〒 599-8570 大阪府堺市中区学園町 1-2
Particle in Cell (PIC) simulation is a very powerful tool to study the sheath structure. But
careless choice of geometry parameter easily lead to unstable and unrealistic sheath behavior. In
this work, 1D simulation results are compared with previous reported kinetic model results. In
order to apply PIC simulation to more complicated probe system, 2D simulation and its stability
are also studied. From these study, importance of careful simulation model for simultaneous
simulation of sheath fine structure and large scale bulk plasma is confirmed.

2. Kinetic theory of the sheath
Here we restrict the situation to stationary
1-dimensional case for simplicity [3]. Plasma
particles is injected from the plasma boundary and, if they reach the right wall boundary, they are absorbed perfectly. Potential atplasma boundary is set to be zero and negatively biased externally at wall boundary.

System length from kinetic model

System length/Debye length

1. Introduction
In order to study the sheath, various boundary condition of not only the solid surface side
but also the plasma upﬂow side must be considered. In analytical model, moreover, calculation geometry must be divided into quasineutral presheath and collisionless sheath, whosesolutions must be connected at a singular point( so-called sheath boundary ) in a rather artiﬁcial manner. Moreover analytical model
often assumes that electrons obey MaxwellBoltzmann relation and that ion’s motion canbe described by so-called free-fall model, which
is not realistic in many plasma ﬁelds.
Particle simulation is a very powerful tool
to study the sheath problem, since it is freefrom the above mathematical procedure. But
careless choice of geometry parameter easily
lead to unstable and unrealistic sheath behavior. In the present study, Particle in Cell
( PIC ) simulation by using Berkeley code (XOOPIC ) is done [1]. Firstly, 1D simulation
results are compared with the kinetic model
results. In order to apply PIC simulation tomore complicated probe system[2], 2D simulation and its stability are also studied.

30

2010No sheath

1.26

1.28

1.3

E–field/(Te/Debye length) at boundaryFig.1: Condition for sheath establishment.
Horizontal axis is the Absolute value of electric ﬁeld (|Es |) at source boundary. Verticalaxis is the length where the potential reaches
assigned solid boundary value.

We also neglect collision and particle sources( sinks ) in the simulation geometry. Distribution function of particle spices j obeys 1D
steady collisionless Boltzmann equation.
vx∂fj (x, vx ) −qj dΦ ∂fj (x, vx )+= 0,∂x
mj dx
∂vx

(1)

By integrating the solution of eq (1) withvx , density is obtained as the function of normalized potential ϕ = e(Φs − Φ)/Te and τ =Ti /Te . Potential ϕ obeys normalized Poissonequation.
−dE
dϕ= E(s),= ρ(s) = (ni −ne )/nes , (2)dsds

where s = x/λD .

3. Simulation result
Firstly, a very simple simulation geometryis examined in Cartesian (x–y) Coordinate tocompare with the kinetic model results. Thesimulation domain is a rectangle, and homogeneity along y-direction is kept. Left boundary at x = 0 is the plasma source and its potential is assumed to be zero. Right boundary
at x = Lx is the conductor, whose potential iskept to be -100[V].
If the domain size (Lx ) is not so large, theXOOPIC simulation reaches the steady stateand S–shape potential proﬁle is obtained. When
Lx increases, however, the situation changes.
Potential proﬁle shows small ﬂuctuations andits amplitude becomes large with Lx . Eventually, for Lx /λD > 30, potential proﬁle oscillates and simulation does not reach the steadystate solution. In this simulation, plasma density is determined by the input value of sourceboundary current. Even keeping Lx constant,by increasing this current , Debye length λD
becomes small and sheath structure becomes
unstable. Fig. 2 shows this tendency.
Recently, very complicated probe heads areproposed to extend the parameters to be measured. The sheath structure around these probehead also becomes complicated. Figure 3 showsthe example of 2D-potential proﬁle around therim region of cylindrical probe head, which is
set at right-bottom of this ﬁgure. Although
the distance between probe and boundary inthe z-direction is small than the critical value
obtained with 1D simulation, the distance inthe r-direction is somehow large. So large potential oscillations are observed and the sheathis unstable. In order to make optimize theprobe structure[2], the gap in r-direction must
be reduced than this model geometry.

current=0.30[A]current=0.40[A]current=0.55[A]current=0.70[A]current=0.70[A]200Potential

Although eq. (2) is very stiﬀ, solutions
with diﬀerent boundary electric ﬁeld strength
Es = ( dϕds ) can be obtained. Potential proﬁleobtained is found to be strongly dependent on
precise value of boundary electric ﬁeld. WhenEs = −1.30, potential reaches ϕw at s ∼ 10.
But when |Es | becomes smaller by only 3%,
ﬂat region of potential grows and simulationgeometry size becomes larger by factor of three.
As shown in Fig. 1, such |Es | corresponds tothe limit of stable sheath formation and, forsmaller |Es | value, no stable sheath potential
can not be obtained.

-20
-40
-60
-80
-100
00.01

0.02

0.03

0.04

0.05

0.06

z[m]

Fig.2: Example of 1D sheath simulation fordiﬀerent plasma density. When source boundary current is large, plasma density is alsolarge.

Fig.3: An example of 2D unstable sheath.
Very large potential oscillation is observed atleft-bottom of the ﬁgure.

Acknowledgements
One of Authors (H.M.) acknowledges Prof. Y.
Tomita of National Institute for Fusion Science and Prof. D.Tskhakaya of Universityof Innsbruck for their advice on sheath kinetic model. The authors thank Plasma Theory and Simulation Group in UC-Berkeley for
their excellent PIC codes

References[1]

J.P.Verboncoeur, A.B.Langdon and N.T.Gladd:Comp. Phys. Comm. 87 (1995) 199.

[2]

H.Matsuura and N.Inagaki: Proc. JSPF annual
meeting, 2010, Sapporo, 01P56 [in Japanese].

[3]

H.Matsuura and Y.Tomita: Proc. SPP24, 2007,Osaka, P2-45.