User Tools


Particle emission characteristics

Particles are emitted from walls for two possible reasons. Either a Particle source is attached to the wall or a reaction involving reaction products takes place.

In both cases, the default behaviour is to emit particles with a Maxwellian velocity according to the given wall temperature and a cos(θ) like angular distribution. A Maxwellian distribution in 3D space is characterized by following distribution function in the velocity $v$,

\large$ f(v)dv = \sqrt{\frac{2 m^3}{\pi k^3T^3}}v^2\exp\left(-\frac{mv^2}{2kT}\right),$

where $T$ is the temperature, $k$ the Boltzmann constant, and $m$ the molecular mass. When sitting at a wall surface, the velocity distribution function has to be multplied with the particle velocity $v_z$ normal to the surface, since the arrival of fast particles is statistically emphasized. In the software implementation, the particle velocities of diffusely reflected or generated particles with Maxwellian distribution are generated as follows in cylindric coordinates $(r, \phi, z)$ with the z coordinate aligned parallel to the surface normal:


  \begin{eqnarray}
  v_z &=& \sqrt{-\frac{2kT}{m}\ln(R\#)\nonumber\\
  v_r &=& \sqrt{-\frac{2kT}{m}\ln(R\#)\nonumber\\
  \phi&=& 2\pi R\#\nonumber\\
  \vec{v}_{\rm Emitted} & = & v_z \hat{n} + v_r\cos(\phi)\hat{t_1} + v_r\sin(\phi)\hat{t_2},\nonumber
  \end{eqnarray}

with R# being a random number between $[0, 1[$ and $\hat{n}, \hat{t_1}, \hat{t_2}$ being the surface normal and normalized tangential vectors. If this algorithm is performed for reflected particles in a closed gas volume, the gas will be everywhere in thermal equilibrium with the wall and will have a Maxwellian distribution according to the wall temperature.

It is possible to modify the behaviour of generated particles by supplying an additional command just after the last wall command involving reaction products (either after the last source or reaction).

Maxwellian distribution with alternative temperature

If required it is possible emitting particles with Maxwellian distribution with a different temperature than specified in the wall properties. For this purpose, there is a function set_emission_maxwell which takes the desired temperature as argument:

set_emission_maxwell("species", temperature);

Source of sputtered particles

A sputter process results from random ballistic collisions between atoms inside the target material. The distribution of sputtered species is non-thermal. A thorough discussion on that topic is given in the book “Sputtering by Particle Bombardment” by Behrisch/Eckstein 1). The distribution function of sputtered particles has been discussed by Thompson2) yielding following relation:

\large$ f(E,\theta) \propto \frac{\cos(\theta)}{E^{2-2m}\left(1+U_b/E\right)^{3-2m}}.

Here, m is a parameter describing the interatomic potential between target atoms. Assumption of the hard-sphere approximation leads to m–>0. The software implementation uses this simplification within the function set_emission_thompson together with an angular distribution function $f_\alpha(\theta)$ depending on a parameter $\alpha$,


\begin{eqnarray}
 f(E,\theta) &\propto& \frac{E}{\left(E+U_b\right)^3}\times f_\alpha(\theta)\nonumber\\
 f_{\alpha}(\theta) & = & \left\{\begin{array}{ll}
   \cos^{\alpha}(\theta); & \alpha \ge 0\\
   \cos(\theta) \left(1+\alpha\cos^2(\theta)\right); & \alpha < 0 \end{array}\right. \nonumber
\end{eqnarray}

The angular distribution for negative values of $\alpha$ is an empirical formula for off-axis sputtering characteristics taken from Yamamura, 19913). It could be relevant for polycrystalline target materials and in case of low ion energies4). For performing ballistic simulation of the sputter process in targets, there is a free software called SRIM, which can produce sputter yields, velocity and angle distribution of sputtered particles as well as implantation profiles for arbitrary target compositions and ion impact. The software implementation includes both formula for the angular distribution function, which are selected according to the parameter $\alpha$.

As can be seen in the following graph 5), the distribution of sputtered particles decays much slower with energy than e.g. a Maxwellian distribution with the same peak energy (The peak energy in case of a sputter distribution is located at $E_{\rm peak} = U_b/2$).

One conceptual problem of the above distribution function by Thompson is that there is no upper limitation regarding the kinetic energy of sputtered species. Generation of particles with arbitrarily high kinetic energy is problematic since they potentially violate numerical constraints with respect to time step and cell spacing.

A practical solution is based on the maximum fraction $\lambda$ of energy which can be transferred from an impinging ion with mass $m_i$ to a target atom with mass $m_t$ by perpendicular collision6),

\large $\lambda = \frac{4m_{i} m_{t}}{\left(m_i+m_t\right)^2}.$

With $E_i$ as the energy of impinging ions it is reasonable to assume that the energy of sputtered atoms can be no higher than $\lambda\times E_i$. Thus, in the DSMC / PICMC software implementation, $\lambda\times E_i$ is used as upper limit of the kinetic energy of sputtered particles. The accordant function for the Thompson distribution is

set_emission_thompson("species", Ub, alpha, E_ion, m_ion);

The parameters are as follows:

  • “species”–> the name of the sputtered species
  • Ub –> The binding energy
  • alpha –> The parameter $\alpha$ determining the angular distribution
  • E_ion –> The kinetic energy of impinging ions
  • m_ion –> The relative molecular mass of impinging ions

The latter two parameters are used to compute the above mentioned recoil factor $\lambda$. The simplified Thomson distribution can be analytically integrated, thus in the practical implementation the energy is taken via

$ E = \frac{U_b}{1-\sqrt{R\#}} - U_b

Again, R# is a random number within[0,1[. The emission angle is created via

$ \theta = \left\{\begin{array}{ll}
  \arccos\left(\left(1-R\#\right)^{\frac{1}{\alpha+1}}\right); & \alpha \ge 0\\
  \arccos\left(\sqrt{\frac{\sqrt{1+\alpha(2+\alpha)(1-R\#)}-1}{\alpha}}\right); & \alpha < 0
  \end{array} $

If the random number yields $E>\lambda\E_i$ then the value is dropped and the calculation is repeated for another random number R#.

Note: In former versions (<2.4.2.8) the Thompson distibution function was generated via

set_emission_sputter("species", Ub, alpha);

Here, the cut-off energy was internally hard-wired as $E_{\rm cutoff} = 100\times U_b$. This function still works but we recommend using the newer function set_emission_thompson instead, where the cut-off energy can be freely specified. Especially in case of sputtered species deposited onto a substrate it turned out that unrealistically high values of the cut-off energy can lead to significant deviations in their mean energy.

Uniform energy distribution within a given interval

For a uniform distribution within an energy interval between E0…E1, the following function can be used:

set_emission_energy("species", E0, E1, alpha);

The angular distribution as a function of the parameter $\alpha$ is the same as in case of sputtered particles (see previous section).

Example: Niobium sputter source

add_source("Nb", 50, "sccm");          # Constant sputter rate equivalent to 50 sccm particle flux
set_emission_sputter("Nb", 7.39, 1.5); # Sputter emission with binding energy 7.39 eV and cosine index 1.5

Particle emission distribution by user defined histogram

If the built-in particle emission profiles are not sufficient, there is also the possibility of implementing a user-defined distribution function via histogram data files. The histogram data files must have the same format as the files resulting from sampling energy and angular distribution, i.e. four columns comprising

  • Energy [eV]
  • Theta [°]
  • Phi [°]
  • Scalar intensity

There are the following four options:

  1. Energy only distribution function
    set_emission_histogram_E("species", alpha, hist_filename);

    Here, only the first and fourth column of the histogram file (energy and intensity) will be taken into account. For the angular distribution function, the parameter alpha has the same meaning as shown above.

  2. Angular only distribution function
    set_emission_histogram_A("species", temperature, hist_filename, ex, ey, ez);

    Here, the first column of the histogram file is ignored; relevant are columns 2-4 (theta and phi angle as well as intensity). The angles are to be given in degree. Additionally, a unit vector ex, ey, ez is to be specified, this vector defines the direction, where the azimuth angle is zero.

  3. Separated energy and angular distribution functions
    set_emission_histogram_EA("species", hist_E_filename, hist_A_filename, ex, ey, ez);

    Here, energy and angular distribution function are both to be specified and treated independently.

  4. Correlated energy and angular distribution function
    set_emission_histogram_Corr("species", hist_filename, ex, ey, ez);

    This function uses only one histogram file, where energy and angular columns are both included, i.e. energy and angular distributions are correlated.

Particle emission distribution by user defined functions

Instead of providing emission distribution functions as histogram files, they can be also applied as analytic formula. There are the following options:

  1. Energy distribution function
    set_emission_function_E("species", alpha, eMax, ne, formula); 

    The energy distribution function for the emitted “species” is defined on an energy interval 0…eMax with ne being the number of slots. The energy unit used to specify the range is eV. The parameter alpha specifies the angular emission profile as specified here. The analytic function is specified via parameter formula, which can contain an analytical formula with E being the energy in eV.

  2. Angular distribution function
    set_emission_function_A("species", T0, nTheta, nPhi, formula, ex, ey, ez);

    This specifies an angular distribution function while assuming an Maxwellian energy distribution function according to the temperature T0. The resolution of the angular distribution function is defined by the number of slots, nTheta for the polar angle and nPhi for the azimuth angle. The analytical function in parameter formula may contain the variables THETA (range: 0-90° in degree), theta (range: $0\ldots\pi$), PHI (range: 0-360° in degree) or phi (range: $0\ldots 2\pi$). The parameters ex, ey, ez represent a vector which defines the direction where $\phi=0$.

  3. Independent energy and angular distribution functions
    set_emission_function_EA("species", eMax, ne, formula_E, nTheta, nPhi, formula_theta_phi, ex, ey, ez);

    This function - as a combination of the two previous functins - allows to set independent energy and angular distributions. The energy distribution function for the emitted “species” is defined on an energy interval 0…eMax with ne being the number of slots. The energy unit used to specify the range is eV, and the analytical formula as a function of E is specified via parameter formula_E. The resolution of the angular distribution function is defined by the number of slots, nTheta for the polar angle and nPhi for the azimuth angle. The analytical function parameter formula_theta_phi may contain the variables THETA (range: 0-90° in degree), theta (range: $0\ldots\pi$), PHI (range: 0-360° in degree) or phi (range: $0\ldots 2\pi$). The parameters ex, ey, ez represent a vector which defines the direction where $\phi=0$.

  4. Correlated energy and angular distribution function
    set_emission_function_Corr("species", eMax, ne, nTheta, nPhi, formula, ex, ey, ez);

    This function allows to specify a correlated energy and angular distribution function. The energy distribution function for the emitted “species” is defined on an energy interval 0…eMax with ne being the number of slots and the energy unit in eV. The resolution of the angular distribution function is defined by the number of slots, nTheta for the polar angle and nPhi for the azimuth angle. The distribution function specified in formula can contain all variables E, THETA, PHI, theta, phi as specified above.

Specifying energy loss functions

In case of e.g. reflected neutrals from ion impact or reflected high energetic electrons, where a certain energy loss occurred in the material, the energy distribution of the reflected species is not given as a function of absolute energy. Instead, there is a distribution function for either the absolute or the relative energy loss. Such cases can be specially treated:

  1. Relative energy loss function
    set_energy_loss_rel("species", formula);

    A relative loss function occurs e.g. for fast reflected neutrals resulting from ion impact on a target. Here, the ratio between impact and emitted energy depends e.g. on a mass ration between ion and target atom and thus, the loss energy distribution is defined on a relative scale E=0…1. This function uses a resolution of 1000 slots of the relative energy scale by default. The parameter formula may depend on the variable E, which - in this case - is not the absolute energy but a relative one (E=0…1), where E=1 corresponds to the impact energy.

  2. Absolute energy loss function
    set_energy_loss_abs("species", eMax, ne, formula);

    This realizes an absolute loss function. One scenario may be a high energetic electron impinging on a surface thereby inducing some Auger excitation or other energy loss mechanism in the material which would substract an absolute loss energy from the electron. In this case, the energy range is E=0…eMax with the unit in eV, and the number of slots is specified in the nE parameter.

In these energy loss functions, the direction of the reflected particle is currently set specularly with respect to the impinging particle (even if e.g. the species changed from ion to neutral). In future versions this might be modified e.g. by specifying an angular cone for the reflected particle.

Example: In plasma simulations, the creation of fast reflected neutrals by ion impact can be realized as follows:

add_plain_reaction("Arplus", 1, 1, ["Ar", "e"], [1, 0.1]);
set_emission_function_E("e", 1.0, 10, 100, exp(-((E-5)^2)/2));
set_energy_loss_rel("Ar", E>0.35);

In this example, an Arplus ion impinging onto a surface results in an electron with a secondary electron yield of 0.1 (10%) and a neutral Ar atom7) with a yield of 100%. From the reflected neutral (Ar) the relative loss energy distribution is set to E>0.35. In the relative scale E=0…1 this evaluates to 1 of E>0.35 and 0 if E⇐0.35. It means that the energy of the reflected Ar is reduced by 35% or more in comparison with the impact energy of the ion. In other words: The energy loss is equally distrbuted between 35% and 100% (The energy distribution of the secondary electron is set as Gaussian function centered around 5 eV, but this is not important, here).

1)
R. Behrisch, W. Eckstein: Sputtering by Particle Bombardment, Springer Science & Business Media, page 135, Table 1, https://books.google.de/books?id=K2Oh8QQ1Kf0C&pg=PA135
2)
M.W. Thompson, Physical Mechanisms of Sputtering, Physics Reports 69 (1981) 335-371 –> see equation (12).
3)
Yamamura, Y.; Takiguchi, T.; Ishida, M.: Energy and angular distributions of sputtered atoms at normal incidence, In: Radiation Effects and Defects in Solids 118 No. 3 (1991) p. 237-261; see equation (13) on page 245.
4)
G.K. Wehner, D. Rosenberg, Angular Distribution of Sputtered Material, Journal of Applied Physics 31/1 (1960) 177-179
5)
the graph is adapted from A. Pflug, M. Siemers, T. Melzig, D. Rademacher, T. Zickenrott, M. Vergöhl, Numerical optimization of baffles for sputtering optical precision filters, Surface and Coatings Technology 241 (2014) 45-49.
6)
M.W. Thompson, The energy spectrum of ejected atoms during the hig energy sputtering of gold, Philosophical Magazine 18 (1968) 377-414.
7)
If the neutral density of Ar is very high with respect to ions and electrons, it might make sense to define a special species such as ArFast in order to distinguish the fast reflected neutrals from the background gas