Transcript Document

Charged particle kinetics by
the Particle in Cell / Monte
Carlo method
Savino Longo
Dipartimento di Chimica dell’Università di Bari and IMIP/CNR
The system under examination
A gas can be ionized under non equilibrium conditions (too low
temperature for equilibrium ionization) with constant energy
dissipation, like in electric discharges, photoionized media,
preshock regions, and so on.
The result is a complex system where the nonlinear plasma
dynamics coexists with chemical kinetics, fluid dynamics,
thermophysics and chemical kinetics issues
Basic phenomenology
The gas is only weakly ionized
Molecules are only partially dissociated and exhibit their chemical
properties
The electron temperature is considerably higher (about 1eV)
than the neutral one (< 1000K)
Velocity and population distributions deviate from the equibrium
laws i.e. Maxwell and Boltzmann respectively
Items to be included in a
comprehensive model
Plasma dynamics
Neutral particles and plasma interaction
Chemical kinetics of excited states
I
Plasma dynamics
The problem of plasma dynamics
The charged particle motion is affected by the electric field, but
the electric field is influenced by the space distribution of
charged particles (space charge)
F  qE / m
2
   4 

Particle in Cell (PiC) method
The method is based on the simulation of an ensemble of mathematical
“particles” with adjustable charge which move like real particles and a
simultaneous grid solution of the field equation
Integration of equations
of motions, moving particles
E field
Grid to particle
Interpolation
Dt
Particle to grid
Interpolation
solve Poisson Equation
for the electric potential
Charge density
Ideal plasma
1
V
g

3
3
nsimD N simD
1
Vlasov equation
eE



v



r
v  f (r , v, t )  0

m
 t

E  
 0 2  ion  e  fd 3v
Particles propagate the initial condition
moving along characteristic lines of
the Vlasov equation
Particle/grid interpolation: linear
q(iDx)   q p S ( x p  iDx)
p
Dx  D
Particle move: “leapfrog”
qE
Dt
m
x  x  v Dt
v  v 

Dt  1/  pl
Plasma oscillation
II
Plasma dynamics
+
Neutral particles and plasma interaction
Vlasov-Boltzmann equation
eE


 v  f (r, v, t )  Cf
  v  r 
m
 t

E  
 0 2  ion  e  fd 3v
Cf  1   /  max  f 
1
 max
3
d
 vpvv f ( v)
pvv   d 3 wd 3 w (v, w; v, w) gF (r, w )
 (r, v )   d 3v pv v
Lagrangian Particles as propagators
Vlasov equation
Initial d moves along
characteristic lines -->
deterministic method
(PIC)
Vlasov/Boltzmann equation
medium
tcoll ?
v’ ?
event
“free flight”
Dispersion of the initial d -->
“choice” -->
stochastic method
(MC)
Statistical sampling of the linear collision operator
(1) Sampling of a collision
partner velocity w from the
distribution F(r,w)/n
 r,v 
1
Af  

1 
f r ,v 

max 
max
pvv 



d 3v pvv f v
d 3 wd 3 w  ( v,w;v, w)| v w| F(r ,w)
(2) rejection of null collisions with
probability 1-ng(g)/ max
(3) kinematic treatment of the
collision event for the
charged+neutral particle system
Test particle Monte Carlo
A ‘virtual’ gas particle is generated as a
candidate collision partner based on the
local gas density and temperature.
The collision is effective with a probability
n gasg
max(n gasg)
For an effective collision the new velocity
of the charged particle is calculated
according to the conservation laws and the

differential cross section
A random time to the next candidate
collision is generated
Preliminary test: H3+ in H2
17
10
1
10
0
15
14
T=600 K
13
T=300 K
mean energy (eV)
K 0 (cm 2V -1 s -1 )
16
T=600 K
10 -1
T=300 K
12
10 -2
11
10 1
10 2
E/n(Td)
reduced mobilities of H3+ ions as
a function of E/n compared with
experimental results of Ellis2 (dots)
2H.
10
100
E/n(Td)
mean energy of H3+ ions as a
function of E/n
W. Ellis, R. Y. Pai, E. W. McDaniel, E. A. Mason and L. A. Vieland, Atomic Data Nucl. Data Tables 17, 177 (1976)
Example: H3+/H2 transport* in
a thermal gradient = 500 K/cm, costant p = 0.31 torr
f(x,y,0)=d(x) d(y) d(x)
E/N=100 Td
* only elastic collisions below about 10eV
7ms no field
7ms with
E field
1,2
7ms
1
f(x,y)
0,8
-1
h(y) (cm )
f(y)
no field
0,6
E field
0,4
0,2
0
-2
-1
0
y (cm)
1
2
Particle in Cell with Monte Carlo Collisions (PiC/MCC) method
Monte Carlo Collisions
Integration of equations
of motions, moving particles
E field
Grid to particle
Interpolation
Dt
Particle to grid
Interpolation
space charge
solve Poisson Equation
for the electric potential
Making the exact MC collision times compatible
with the PIC timestep
After R.W.Hockney, J.W.Eastwood, Computer Simulation using Particles, IOP 1988
Plasma turbulence due to charge exchange in Ar+/Ar
(collaboration with H.Pecseli , S. Børve and J.Trulsen, Oslo)
2 component (e,Ar+) 1.5D PIC/MC
106 superparticles
vx
t=0
Initial beam:
 = 4 1013 m-3
< > = 1eV
T = 100 K
L = 0.05 m
Ar background:
T = 100K, p= 0.3torr
The electron density is calculated as a
Boltzmann distribution, this produces a
nonlinear Poisson equation solved iteratively
x
1  2 
4 e
nion  ne0 exp(e / kTe )
vx
electrostatic
repulsion
inertia
collisions
x
The collisional production of the second (rest) ion beam can lead to a
two stream instability
Two stream instability
The propagation of two charged particle beams in opposite
directions is unstable under density/velocity perturbations and can
lead to plasma turbulence
v
r
(v2  v1 )   pl
vx (m/s)
1
10 (log)
Simulation time. 2 10-5 s
4000
2000
0
-2000
0.050
0.100
0.150
x (m)
0.200
0.250
Capacitive coupled, parallel plate
radio frequency (RF) discharge
V RF sin( 2 RF t )

negative
charge
strong oscillating field
regions (sheaths)
electrons
electron density
ambipolar potential
energy well = -e
negative
charge
Simplified code implementation for nitrogen
2 particle species in the plasma phase: e, N2+
more than one charged species
Selection of the collision process based on the cross section database
Process probability = relative contribution to the collision frequency
Particle position/energy plot
)
10
-3
electron and ion density (m
V = 500 V, p = 0.1 torr, f =13.5 MHz
16 rf
10
15
ions
10 14
electrons
10
13
0
0.01
0.02
position (m)
0.03
0.04
electron/ion mean energy (eV)
Vrf = 500 V, p = 0.1 torr, f =13.5 MHz
10
electrons
1
ions
0.1
0
0.01
0.02
0.03
position (m)
0.04
electron/ion drift velocity (m/s)
Vrf = 500 V, p = 0.1 torr, f =13.5 MHz
4
1 10
5000
ions
0
-5000
electrons
4
-1 10
0
0.01
0.02
0.03
position (m)
0.04
III
Plasma dynamics
+
Neutral particles and plasma interaction
+
Chemical kinetics of excited states
Kinetics of excited states
e A e A *
A* A  h
A *B  A (h)  B
e  A*  e(h)  A
e  A *  2e  A 
A 2* A 2 (v)  h
A *B  A  B
e  A 2 (v)  A   B

Numerical treatment of state-to-state chemical
kinetics of neutral particles (steady state)
(1) gas phase reactions:



c1
Nc
 rc X c  
Nc
c1

1rN
rc X c 
E.g.: H  H 2 (X ,v )  H  H 2 (X , v )
r
 nc x
2
are included by solving:
Dc

x2

r rc  rc kr fe t c nc
rc
(2) gas/surface reactions:
1
E.g.: H( wall)  H 2
2
A s   r1A1   r2 A 2  ...

are included by setting appropriate boundary conditions


Ds s  

r
 rss 
rs

 rs rss
s 
1 8KT
s
4 m s
Boundary
Conditions
surface
reactions
(wall)
Poisson Equation
absorption,
sec.emission
Reaction/Diffusion
Equations
electric
field
eedf
electr./ion
density
gas composition
space
charge
Charged Particle
Kinetics
Monte Carlo Collisions
Chemical
kinetics
equations
Integration of equations
of motions, moving particles
E field
Grid to particle
Interpolation
Particle to grid
Interpolation
solve Poisson Equation
for the electric potential
Space
charge
code implementation for hydrogen
5 particle species in the plasma phase: e, H3+, H2+, H+, H16 neutral components: H2(v=0 to 14) and H atoms
Charged/neutral particle collision processes
electron/molecule and electron/atom elastic, vibrational and electronic
inelastic collisions, ionization, molecule dissociation, attachment, positive
ion/molecule elastic and charge exchange collisions, positive elementary
ion conversion reactions, negative ion elastic scattering, detachment, ion
neutralization
Schematics of the state-to-state chemistry for neutrals
e + H2(v=0)  e +H2(v=1,…,5)
e + H2(v=1,…,5)  e +H2(v=0)
H2(v) + H2(w)  H2(v-1) + H2(w+1)
H2(v) + H2  H2(v-1) + H2
H2(v) + H2  H2(v+1) + H2
H2(v) + H  H2(w) + H
e + H2(v=0,…,14)  H + He + H2(v)  e + H2(v’) (via b1u+, c1u)
e + H2  e +2H(via b3u+, c3u, a3g+, e3u+)
e + H2  H + H+ + 2e
e + H2  H2+ + 2e
H2+ + H2  H3+ + H (fast)
H2(v>0) – wall  H2(v=0)
H – wall  1/2 H2(v)
e + H  2e + H+
e + H2  e + H + H(n=2-3)
e + H-  2e + H
secondary ions from:


H2  H2  H 3  H
charged particle density
1015
Simulation parameters:
Tg = 300 K
Vrf = 200 V
p = 13.29 Pa (0.1 torr)
rf = 13.56 MHz
L = 0.06 m, Vbias = 0 V
v = 0.65, H = 0.02
number density (m -3)

H
1014
-
H3 +
eH+
1013
H2 +
1012
0
0,01
0,02
0,03
position (m)
primary positive ions
0,04
0,05
0,06
relatively low
T01 (~1000K)
1022
number density (m -3)
1020
0.6
1.2
1.8
2.4
1018
cm
cm
cm
cm
1016
1014
1012
0
2
4
6
8
10
12
14
vibrational quantum number
plateau due to
radiative EV processes
eedf
100
0.6
1.2
1.8
2.4
3.0
10-1
eedf (eV -3/2)
10-2
cm
cm
cm
cm
cm
10-3
10-4
10-5
10-6
0
5
10
15
20
25
energy (eV)
30
35
40
Double layer
O. Leroy, P. Stratil, J. Perrin, J. Jolly and P. Belenguer, “Spatiotemporal analysis of the double layer formation in
hydrogen radio frequencies discharges”, J. Phys. D: Appl. Phys. 28 (1995) 500-507
Bias voltage
p = 0.3 torr
L = 0.03 m
H = 0.0033
V = 0.02
A. Salabas, L. Marques, J. Jolly, G. Gousset, L.L.Alves, “Systematic characterization of low-pressure
capacitively coupled hydrogen discharges”, J. Appl. Phys. 95 4605-4620 (2004)
Conclusion
A very detailed view of the charged particle kinetics in
weakly ionized gases can be obtained by Particle in
Cell simulations including Monte Carlo collision of
charged particle and neutral particles.
Items to study in the next future
(students)
Charge particle kinetics in complex flowfields
Collective plasma dynamics in shock waves
Development of new MC methods for electrons
matching the time scale for electron heating
….