Title: Computational Plasma Physics
1Computational Plasma Physics Kinetic modelling
Part 1 2 W.J. Goedheer FOM-Instituut voor
Plasmafysica Nieuwegein, www.rijnh.nl
2What are kinetic methods and when do they apply
Kinetic methods retain information on the
velocity distribution (hydrodynamic/fluid methods
first integrate over velocity space) Needed when
distribution is non-Maxwellian Kinetic methods
are to be preferred when ?mfpgt L or
?collgtT ?mfp and ?coll depend on densities and
cross-sections But what are L and/or
T?? Examples from (plasma) physics?
3Kinetic models non-Maxwellian
Collisions electrons mainly with neutral
species Low degree of ionization Effective
cooling of parts of the energy distribution
function Counteracted by Coulomb collisions at
high degree of ionization (gt10)
E
4Variations in space and time
Boundary layers Transition layers
Transient phenomena Switching on Modulation
5Power modulated discharges
Modulate RF voltage (50MHz) with square wave (1 -
400 kHz)
Observation in experiments UU) optimum in
deposition rate
RF
6Modulated discharge in SiH4
Results from a PIC/MC calculation Cooling and
high energy tail
7Examples of Ion Energy DF at grounded electrode
From Th. Bisschops, Thesis TU/e, 1987
Interaction between E-field and ion motion does
not result in a shifted Maxwellian
8Kinetic models strong spatial variation
Very low pressures L size of vessel (applies
for e,i,n) Space charge boundary layer L
Debye length (applies for e,i) Micro-structures
(etched trenches) L size of structure
(applies for e,i,n) Shocks L extension of
shock (applies for e,i,n)
There may be a difference between momentum loss
and energy loss
9Kinetic models strong temporal variation
Microwave discharges / high frequency RF
discharges (applies for e,(i)) Start-up of
discharges (applies for e,i)
There may again be a difference between momentum
loss and energy loss
10Methods based on direct solution of the Boltzmann
equation
Tricks to solve BE Use symmetry if present
Expand f in some small parameter
11A method especially suitable for electrons
Electrons have a low mass ? high momentum loss
in collisions ? energy loss in inelastic
collisions Elastic scattering ? redistribution
over a sphere in velocity space Small deviation
from isotropic f in the direction of the average
velocity Therefore expansion in Legendre
polynomials Pn(cos ?) with ? the angle between
average and actual velocity f f0(v)
f1(v)P1(cos ?) f2(v)P2(cos ?). Note Ampli
tudes depend on absolute value velocity and
vary in space and time
12An example f0f1cos(?)
f0f1cos(?)
13How to calculate the amplitudes fn
Elastic collisons move electrons from outside in
Transport from Neighbouring volume
V
Vdv
Shift in velocity due to electric field (eEDt/m)
Net change in Dt
14How to calculate the amplitudes fn
Balance of number of particles in shell between v
and vdv in dxdydz Transport in real
space Transport in velocity space Effect of
collisions
15How to calculate the amplitudes fn
Balance of momentum in shell between v and vdv
in dxdydz Transport in real space Transport in
velocity space Effect of collisions
Note that f1 is a vector, directed along the
average velocity (?0)
16Cutting off at f1 The Lorentz approximation
For elastic collisions with atoms/molecules, with
mass M
17Special case homogeneous, steady state
Temperature gas is zero Constant electric field,
average velocity (f1) along E
18Special case homogeneous, steady state
Solution Backward integration, tri-diagonal
system
19Special solutions
Reduced electric field
r-1 s2 Maxwell r 0 s4
Druyvesteyn Druyvesteyn has less
energetic electrons
20Inelastic collisions
- Couple parts of the distribution function that
are far apart - Example Excitation
- electron looses excitation energy (a few to gt10
eV) - electron is set back in velocity
Source proportional to vf0(v)NM?exc(v) Same
holds for ionization Energy new electrons to be
specified
21Inelastic collisions two T distribution
Noble gases have high first excitation energy For
lower energies only elastic energy losses slow
decay of f with v For higher energies large
energy losses fast decay of f with v Resulting
distribution is characterised by two
temperatures
Eexc
Eion
22Inelastic collisions ionization
In ionization Eion is lost Suppose remaining
energy equally divided How many electrons arrive
between Ui-du and Uidu
- Ui-du lt (U-Eion)/2 lt Uidu
- Eion2Ui-2du ltU ltEion2Ui2du
- So factor 2 from energy range factor 2 from new
electron - 4f0(u)u1/2?ion(u)
In steady state problems new electrons
neglected, Usually this has only a minor
influence
23An example, SiH4/H2, with inelastic collisions
EEDFs with 4eV av. Energy in SiH4/H2
non-Maxwellian
24Some quantities (assuming f0 normalized)
25Use of this approach in modelling
Local field approximation Everything expressed
in local E/N-field mobility and diffusion
coefficients reaction rates (ionization,
excitation) average energy Mean energy
approximation Use solution for various
E/N-fields to construct table (mobility,
diffusion coefficient, rates) all as a function
of the average energy (cf. table as function of
temperature for Maxwellian f) Use fluid energy
balance to obtain av.energy in simulation
26Use of the mean energy approximation
Homogeneous gas of given composition,
Nb1...bn EEDF from Boltz.Eqn. Homogeneous
electric field, constant in time
Mobility (?e), Diffusion (De) EEDF Average
energy (? ? 1.5 kTe) Reaction rates for
processes Kproc (E,Nb1,Nb2,..Nbn)
Combine results in table for Kproc (?) , ?e(?),
De(?)
27Modelling the electrons
Look-up table
28One step further time dependent E-field
Important characteristic times Loss of
momentum goes very fast f1 is in
equilibrium with E-field Loss of energy Only
fast in case of inelastic collision
f0 can be out of equilibrium Example reaction
of f0 in SiH4/H2 E0cos(?t) behavior depends
on ratio ? and loss frequencies
29High frequency smaller excursion f0
Collision frequency ? pressure Therefore
normalized to 1 Torr 10 SiH4, 90
H2 EEmcos(?t), f0 at Em, Em/2, 0, -Em/2 (1,2,3,4)
Energy loss
Momentum loss
From Capitelli et al. Pl. Chem. Plasma Proc. 8
(1988) 399-424
30Time dependent, spatially inhomogeneous E field
Is possible in principle, but More than 1
spatial dimension would take too much CPU
time Really steep gradients (sheaths) require fn
with ngt1 Solution Monte Carlo methods Account
in principle for all effects
31Example vf0 in Nitrogen
E3.6104(x/L)5(1?0.8sin(?t)), -LltxltL ?2?80 MHz
32Electron-electron collisions
Electrons efficiently share energy in elastic
collisions ? Collisions try to establish
Maxwell distribution
More sophisticated operators conserve momentum
and energy
33Monte Carlo methods
Principle Follow particles by - solving
Newtons equation of motion - including the
effect of collisions - collision an event
that instantaneously changes the velocity Note
The details of a collision are not modeled
Only the differential cross section effect
on energy is used Examples Electrons in a
homogeneous electric field Follow sufficient
electrons for a sufficient time Obtain
distribution over velocities etc. ? f0,f1
Positive ions in plasma boundary layer (ions
have trouble loosing momentum)
34Monte Carlo methods Equation of motion
Leap-frog scheme
Alternative Verlet scheme
35Monte Carlo methods B-field
Problem with Lorentz force contains velocity,
needed at time t Solution take average
The new velocity at the right hand side can be
eliminated by taking the cross product of the
equation with the vector
36Monte Carlo methods Boris for B-field
Equivalent scheme (J.P.Boris), (proof
substitution)
37Monte Carlo methods Collisions
Number of collisions NM?tot 1/? per meter. ?
?(x) ?(0)exp(- NM?x) ?(0)exp(-x/?) dP(x)fr
action colliding in (x,xdx)exp(-x/?)(1-exp(-dx/?
))(dx/?)exp(-x/?) ? P(x)(1-exp
(-x/?)) Distance to next collision
Lcoll-?ln(1-Rn) (Rn is random
number,0ltRnlt1) Number of collisions NM?tot v
1/? per second. Time to next collision
Tcoll-? ln(1-Rn)
38Monte Carlo methods Collisions
- Another approach is to work with the chance
- to have a collision on v?t Pcv?t/?
- Ensure that v?tltlt? to have no more than one
collision per timestep - Effect of collision just after advancing
position or velocity - introduces only small error
- When there is a collision
- Determine which one new random number
39Monte Carlo methods Null Collision
Problem Mean free path is function of velocity
Velocity changes over one mean
free path Solution Add so-called
null-collision to make v?tot independent of v
Null-collision does nothing with
velocity Mean free path thus based on Max
(v?tot) Is rather time-consuming when v?tot
peaks strongly
40Monte Carlo methods Null Collision
41Monte Carlo methods Effect of collision
Determine effect on velocity vector Retain
velocity of centre of gravity Select by random
numbers two angles of rotation for relative
velocity Subtract energy loss from relative
energy Redistribute relative velocity over
collision partners Add velocity centre of gravity
42Monte Carlo methods Effect of collision
v1,v2 velocities in lab-frame prior to collision,
w1,w2 in center of mass system
43Monte Carlo methods Effect of collision
A collision changes the size of the relative
velocity if it is inelastic
- A collision rotates the relative velocity
- Two angles of rotation ? ? ????? and ? ? ??????
- usually has an isotropic distribution ?Rn??
- has a non-isotropic distribution
- Hard spheres
44Monte Carlo methods Rotating the relative
velocity
Step 1 construct a base of three unit
vectors Step 2 draw the two angles Step 3
construct new relative velocity Step 4
construct new velocities in center of mass
frame Step 5 add center of mass velocity
45Monte Carlo methods Applicability
- Examples where MC models can be used are
- motion of electrons in a given electric field in
a gas (mixture), see practicum - motion of positive ions through a RF sheath
(given E(r,t))
46Monte Carlo methods Applicability
- Main deficiency of Monte Carlo not
selfconsistent - electric field depends on generated net electric
charge distribution - current density depends on average velocities
- following all electrons/ions is impossible
- Way out Particle-In-Cell plus Monte Carlo
approach
47Particle-In-Cell plus Monte Carlo the basics
- Interactions between particle and background gas
are dealt with only in collisions - this means that PIC/MC is not! Molecular Dynamics
- each particle followed in MC represents many
others superparticle - Note each superparticle behaves as a single
electron/ion - Electric fields/currents are computed from the
superparticle densities/velocities - -But charge density is interpolated to a grid,
so no delta functions
48Particle-In-Cell plus Monte Carlo Bi-linear
interpolation
xs
zi1(i1)?z
xs, qseNs
zs
zii?z
xii?x
xi1(i1)?x
?i?i(xi1-xs)qs/?x ?i1?i1(xs-xi)qs/?x
xjj?x
xj1(j1)?x
?ij?ij(zi1-zs) (xj1-xs) qs/(?x ?z)
49Particle-In-Cell plus Monte Carlo Solution of
Poisson equation
?2
Boundary conditions on electrodes, symmetry, etc.
Electric field needed for acceleration of
particle (bi)linear interpolation, field known
in between grid points
50Particle-In-Cell plus Monte Carlo Full cycle,
one time step
51Particle-In-Cell plus Monte Carlo Problems
Main source of problems Statistical
fluctuations Fluctuations in charge
distribution fluctuations in E average is zero
but average E2 is not ? numerical heating Sheath
regions contains only few electrons Tail of
energy distribution contains only few
electrons large fluctuations in ionization rate
can occur
52Particle-In-Cell plus Monte Carlo Problems
- Solutions
- Take more particles (NB error as N-1/2 ) ,
parallel processing! - Average over a long time
- Split superparticles in smaller particles when
needed - requires a lot of bookkeeping, different weights!
53Particle-In-Cell plus Monte Carlo Stability
Plasmas have a natural frequency for charge
fluctuations The (angular) Plasma
Frequency And a natural length for shielding of
charges The Debye Length Stability of PIC/MC
requires
54Power modulated discharges
Modulate RF voltage (50MHz) with square wave (1 -
400 kHz)
Observation in experiments UU) optimum in
deposition rate
55Modulated discharges
Results from a PIC/MC calculation Cooling and
high energy tail
561-D Particle-In-Cell plus Monte Carlo
Simulation of a dusty argon plasma
Capture cross section
Crystal (2?1010 m-3) 7.5 ?m radius
Scattering Coulomb, truncated at ?d
Void
L/8
L/4
RF
w is energy electron/ion
57Charging of the dust upon capture of ion/electron
The total charge is monitored on the
gridpoints Charge of superparticle is added to
nearest gridpoints Division according to linear
interpolation Local dust charge is total charge
divided by nr. of dust particles This number is
densitydz?a2 For Monte Carlo the maximum ?v is
computed Null-collision is used Note the
difference with the collisions with the
uniform background gas here we construct a
grid-related probability of an event
58Simulation for Argon, 50MHz, 100mTorr, 70V, L3cm
dustfree
with dust
Vd?6V
59Simulation for Argon, 50MHz, 100mTorr, 70V, L3cm
dustfree
with dust
60Simulation for Argon, 50MHz, 100mTorr, 70V, L3cm
Generation of internal space charge layers
An internal sheath is formed inside the
crystal Ions are accelerated before they enter
the crystal This has consequences for the
charging shielding
61Particle-In-Cell plus Monte Carlo What if
superparticles collide?
Example recombination between positive and
negative ions Procedure number of
recombinations in ?t NN-Krec ?t
corresponds to removal of corresponding
superparticles randomly remove negative ion
and nearest positive ion but be careful if
distribution is not homogeneous Again a
grid-based probability can be used A more
sophisticated approach Direct Simulation Monte
Carlo
62DSMC Basics
Divide the geometry in cells Each cell should
contain enough testparticles Newtons equation
as before, but keep track of cell
number Collisions choose pairs (in same cell!)
and make them collide Essential the velocity
distribution function is sum of ?-functions
63DSMC Choosing the pairs
Add null collision Chance of collision of
particle i with j is Pc(Npp/Vcell)Max(v?)?t Ave
rage number of colliding pairs n(n-1)
Pc/2 Select randomly n(n-1) Pc/2 particle pairs
(make sure no double selection) See if there is
no null collision, again with random number,
comparing the real chance for this collision
(vr?) with the maximum Max(v?) Perform the
collision
64DSMC An example
Relaxation of a mono-energetic distribution to
equilibrium 20000 particles, hard sphere
collisions, one cell contains all particles
65Fast when possible, kinetic when needed Hybrid
models
Examples of hybrid models Hydrodynamic ions and
cold electrons Monte Carlo for fast
electrons (tail EEDF) Boltzmann electrons,
Monte Carlo for ions MHD model for plasma,
Monte Carlo neutrals Problems are due to
coupling iterations needed
66B2-EIRENE for Magnum-psi
675 slm H2 Th 2eV Te12 eV 1024 m-2s-1
T profile at inlet
68Recycling at the target