5. Example of application: Monte Carlo simulation of a gas mixture in the PECVD

In this section, we present an example of PP-MCS of collisions and reactions in gas phase of SiH4/H2 mixture used in PECVD process. Some paragraphs have been treated in previous works [21, 24].

#### 5.1 Description of the physical phenomenon

We use a MCS to study collisions and chemical reactions in gas phase of SiH4/H2 mixture used in the PECVD process. In this phase, important reactions have been identified that contribute to the production and the consumption of hydrogen (H), silylene (SiH2), and silyl (SiH3). The hydrogen consumption reactions SiH4 + H ! SiH3 + H2 and SiH3 + H ! SiH2 + H2 are found to play a central role in deciding the distribution of hydrogen [39].The plasma chemistry indicates that H atoms and SiH3 radicals play an important role in the a-Si:H deposition process [40]. Experimentally, it is generally accepted that SiH3 radicals dominate a-Si:H and μc-Si film growth from SiH4 plasmas in the PECVD; it is the key precursor of a-Si:H deposition [41]. The proposed MCS allowed to get the ratio SiH2/SiH3 and mean

How to Use the Monte Carlo Simulation Technique? Application: A Study of the Gas Phase… DOI: http://dx.doi.org/10.5772/intechopen.88559

value of densities of species. It provides information on SiH4 dissociation and on the production of SiH3, H, SiH2, and Si2H6 and other important parameters.

The plasma in the PECVD reactor is weakly ionized. For our study, the mixture gas contains 22% of SiH4 and 78% of H2; the pressure is 100 mtorr, the temperature of the gas ranges from 373 to 723 K, the electron temperature is about 2.5 eV, and the electron density is 3. 108 cm�<sup>3</sup> . The process is considered to be stationary. We take into account electrons and eight neutral species (SiH4, SiH3, SiH2, H, H2, Si2H6, Si2H5, SiH). Reactions taken into account include seven electron-neutral and 14 neutral-neutral reactions. Table 1 shows the 21 reactions and rate constants Kreac. At low temperature, the neutrals interact occasionally with each other and move under the effect of thermal agitation; their velocity distribution function is Maxwell-Boltzmann distribution. Electrons have the mean velocity with kinetic energy Te.

Let Kcons reac and Kprod reac be the rate constants of the consumption and the production of species A. The chemical reaction for the consumption of A is as:

$$a \cdot A + b \cdot B \xrightarrow{K\_{\text{max}}^{\text{corr}}} c \cdot C + d \cdot D$$

And chemical reaction for the production of A is as:

$$a' \cdot A' + b' \cdot B' \xrightarrow{K\_{\text{max}}^{\text{prod}}} c' \cdot A + d' \cdot D'$$


Table 1.

List of gas phase reactions and corresponding rate constants [24].

module EM calculates the electromagnetic fields by solving Maxwell equations. These fields are used as inputs in the module MCS, where the electron density, electron temperature, electron energy distribution function, and electron impact reaction rates can be computed with a Monte Carlo procedure. Subsequently, the module fluid calculates densities and fluxes of the various plasma species (i.e., heavy particles and electrons) with continuity equations and the electrostatic field with Poisson's equation. This electrostatic field is used as input again in the EM. This cycle is iterated until convergence. The schematic of the hybrid model is given in

Theory, Application, and Implementation of Monte Carlo Method in Science and Technology

Schematic of a hybrid model of three modules used to study gas mixtures in the PECVD [33].

To solve statistical physics problems with evolutions as a function of time, kinetic models of MCS (kMCS) are used. Using kMCS, Battaile and Srolovitz [17] described kinetic phenomena of the diffusive motion of a single interstitial atom in a close-packed metal crystal. The motion of the interstitial atom is usually limited to two types: vibration of the atom around the center of the interstitial hole in which it resides and hops to nearest-neighbor interstitial sites. The atom can hop into any of the nearest-neighbor interstitial sites; it executes a random walk. In an MC simulation of this diffusion process, the new position of the interstitial atom is chosen at

Other CVD and PECVD works on MCS are presented in Ref.s [15, 34–38]. They show how MCS methods can study properties of gas mixtures and properties of the

5. Example of application: Monte Carlo simulation of a gas mixture in

In this section, we present an example of PP-MCS of collisions and reactions in gas phase of SiH4/H2 mixture used in PECVD process. Some paragraphs have been

We use a MCS to study collisions and chemical reactions in gas phase of SiH4/H2 mixture used in the PECVD process. In this phase, important reactions have been identified that contribute to the production and the consumption of hydrogen (H),

SiH4 + H ! SiH3 + H2 and SiH3 + H ! SiH2 + H2 are found to play a central role in deciding the distribution of hydrogen [39].The plasma chemistry indicates that H atoms and SiH3 radicals play an important role in the a-Si:H deposition process [40]. Experimentally, it is generally accepted that SiH3 radicals dominate a-Si:H and μc-Si film growth from SiH4 plasmas in the PECVD; it is the key precursor of a-Si:H deposition [41]. The proposed MCS allowed to get the ratio SiH2/SiH3 and mean

silylene (SiH2), and silyl (SiH3). The hydrogen consumption reactions

random from a list of the adjacent interstitial sites.

Figure 4.

Figure 4.

growth of thin films.

the PECVD

126

treated in previous works [21, 24].

5.1 Description of the physical phenomenon

Rate production and consumption for any species A are taken as:

$$R\_A = + \sum\_{K\_{\text{reac}}^{\text{cons}}} K\_{\text{reac}(A', \mathbf{B}')}^{\text{cons}} n\_{A'}^{a'} n\_{B'}^{b'} - \sum\_{K\_{\text{reac}}^{\text{cons}}} K\_{\text{reac}(A, \mathbf{B})}^{\text{cons}} n\_A^a n\_B^b \tag{19}$$

5.2.2 Treatment of elastic and inelastic collisions

target particle j the average collision frequency νij as:

where <sij> is the cross section of the particle j.

The mean free path <λι> of species i is:

The time between two collisions τij is then:

General rules of collision theory are applied:

of energy and momentum for elastic collisions.

• Conservation of total energy as isolated system.

notions of collisions and corresponding parameters [45–48].

bility Prcol,j (nonuniform discrete distribution) given by:

chemical reaction) [21].

129

velocity between the two species i and j.

DOI: http://dx.doi.org/10.5772/intechopen.88559

Let ni and nj be the densities of species i and j in the gas and Vij the relative

How to Use the Monte Carlo Simulation Technique? Application: A Study of the Gas Phase…

<sup>&</sup>lt; <sup>λ</sup>i<sup>&</sup>gt; <sup>¼</sup> <sup>1</sup>

<sup>τ</sup>ij <sup>¼</sup> <sup>&</sup>lt; <sup>λ</sup>i<sup>&</sup>gt; Vij

cies i and j giving products i' and j', the rate constant reaction verifies [45]:

For chemical effective reactions (inelastic collisions) between two reactive spe-

• The new velocities of the colliding particles are calculated using conservation

• Movement of the center of mass and relative motion around the center of mass.

The reader can refer to some fundamental physics books that deal with general

The plasma in the PECVD reactor is weakly ionized. At low temperature, particles interact occasionally with each other and move under the effect of thermal agitation. In reality, only a small fraction of collisions are effective (result in a

In our MCS, after traveling a random walk given by a Gaussian distribution, the

where νij is the neutral-neutral or electron-neutral collisional frequency. The colli-

first chosen particle collides with a second particle (molecule, atom, radical, or electron). The last particle j is randomly chosen according to a (i-j) collision proba-

> Prcol,<sup>j</sup> <sup>¼</sup> <sup>ν</sup>ij P<sup>9</sup> <sup>k</sup>¼<sup>1</sup> <sup>ν</sup>ik

sion theory indicates that the collision between molecules can provide the energy needed to break the necessary bonds so that new bonds can be formed [49]. Particles must have sufficient energy to initiate the reaction (activation energy), so the two chosen particles must have kinetic energy equal to or greater than the barrier energy (Ea) of a

¼ 1 υij

According to the kinetic theory of gases, we have for an incident particle i on a

υij ¼ Vijnj <σij> (21)

kij <sup>¼</sup> <sup>&</sup>lt;σij Vij � � � Vij<sup>&</sup>gt; (24)

nj <sup>&</sup>lt;σij<sup>&</sup>gt; (22)

(23)

(25)

#### 5.2 Description of Monte Carlo simulation technique

#### 5.2.1 Simulation cell and phase coordinates

The MCS is based on binary collisions at the microscopic level. Elastic collisions are between all particles, and inelastic collisions (or effective collisions) are those that result in a chemical reaction. A chemical reaction needs a collision involving at least two particles (atoms, ions, electrons, or molecules). According to kinetic theory, gases consist of particles in random motion. These particles are uniformly distributed in a cell which has a parallelepiped form of sizes Lx, Ly, and Lz (Figure 5). These particles move in a straight line until they collide with other particles or the walls of their container. Dimensions and volume of Monte Carlo cell must take into consideration the mean free path of species.

Let ni be the density of neutral spice i (i = 1,…, 8). The first particle i is randomly chosen according to a probability of neutral species Prsp,I (nonuniform discrete distribution) given by:

$$\text{Pr}\_{\text{sp, i}} = \frac{n\_i}{\sum\_{j=1}^{8} n\_j} \tag{20}$$

The chosen particle takes randomly three components of space in cell ri(xi, yi, zi) according to the normal distribution (nonuniform distribution). It takes also randomly three components of velocity v<sup>i</sup> (vxi, vyi, vzi) according to Maxwell-Boltzmann distribution.

Figure 5. Form of the simulation cell.

How to Use the Monte Carlo Simulation Technique? Application: A Study of the Gas Phase… DOI: http://dx.doi.org/10.5772/intechopen.88559

#### 5.2.2 Treatment of elastic and inelastic collisions

Rate production and consumption for any species A are taken as:

Theory, Application, and Implementation of Monte Carlo Method in Science and Technology

distributed in a cell which has a parallelepiped form of sizes Lx, Ly, and Lz

A0nb<sup>0</sup>

The MCS is based on binary collisions at the microscopic level. Elastic collisions are between all particles, and inelastic collisions (or effective collisions) are those that result in a chemical reaction. A chemical reaction needs a collision involving at least two particles (atoms, ions, electrons, or molecules). According to kinetic theory, gases consist of particles in random motion. These particles are uniformly

(Figure 5). These particles move in a straight line until they collide with other particles or the walls of their container. Dimensions and volume of Monte Carlo cell must take

> Prsp,<sup>i</sup> <sup>¼</sup> ni P<sup>8</sup> <sup>j</sup>¼<sup>1</sup> nj

according to the normal distribution (nonuniform distribution). It takes also randomly three components of velocity v<sup>i</sup> (vxi, vyi, vzi) according to

Let ni be the density of neutral spice i (i = 1,…, 8). The first particle i is randomly chosen according to a probability of neutral species Prsp,I (nonuniform discrete dis-

The chosen particle takes randomly three components of space in cell ri(xi, yi, zi)

<sup>B</sup><sup>0</sup> �<sup>X</sup> Kcons reac

Kcons reac Að Þ ; <sup>B</sup> na

Anb

<sup>B</sup> (19)

(20)

Kcons reac A0 ; <sup>B</sup><sup>0</sup> ð Þna<sup>0</sup>

RA ¼ þ<sup>X</sup>

into consideration the mean free path of species.

tribution) given by:

Figure 5.

128

Form of the simulation cell.

Maxwell-Boltzmann distribution.

5.2.1 Simulation cell and phase coordinates

Kprod reac

5.2 Description of Monte Carlo simulation technique

Let ni and nj be the densities of species i and j in the gas and Vij the relative velocity between the two species i and j.

According to the kinetic theory of gases, we have for an incident particle i on a target particle j the average collision frequency νij as:

$$
\rho\_{\vec{\text{ij}}} = V\_{\vec{\text{ij}}} n\_{\vec{\text{j}}} < \sigma\_{\vec{\text{ij}}} > \tag{21}
$$

where <sij> is the cross section of the particle j. The mean free path <λι> of species i is:

$$<\lambda\_i> = \frac{1}{n\_j < \sigma\_{i\bar{j}}>} \tag{22}$$

The time between two collisions τij is then:

$$
\pi\_{\vec{\eta}} = \frac{<\lambda\_{\vec{\imath}}>}{V\_{\vec{\imath}\vec{\jmath}}} = \frac{1}{\nu\_{\vec{\imath}\vec{\jmath}}} \tag{23}
$$

For chemical effective reactions (inelastic collisions) between two reactive species i and j giving products i' and j', the rate constant reaction verifies [45]:

$$k\_{\vec{\eta}} = \prec \sigma\_{\vec{\eta}}(V\_{\vec{\eta}}) \cdot V\_{\vec{\eta}} > \tag{24}$$

General rules of collision theory are applied:


The reader can refer to some fundamental physics books that deal with general notions of collisions and corresponding parameters [45–48].

The plasma in the PECVD reactor is weakly ionized. At low temperature, particles interact occasionally with each other and move under the effect of thermal agitation. In reality, only a small fraction of collisions are effective (result in a chemical reaction) [21].

In our MCS, after traveling a random walk given by a Gaussian distribution, the first chosen particle collides with a second particle (molecule, atom, radical, or electron). The last particle j is randomly chosen according to a (i-j) collision probability Prcol,j (nonuniform discrete distribution) given by:

$$\text{Pr}\_{\text{col},\dagger} = \frac{\nu\_{\text{ij}}}{\sum\_{k=1}^{9} \nu\_{ik}} \tag{25}$$

where νij is the neutral-neutral or electron-neutral collisional frequency. The collision theory indicates that the collision between molecules can provide the energy needed to break the necessary bonds so that new bonds can be formed [49]. Particles must have sufficient energy to initiate the reaction (activation energy), so the two chosen particles must have kinetic energy equal to or greater than the barrier energy (Ea) of a gas phase reaction. The difference between the kinetic energy of the two particles and the activation energy define the kind of collision (effective or not effective).

The activation energy is given by:

$$E\_a = -k\_B T \ln \left( K\_{\text{reac}} / \nu\_{\vec{\eta}} \right) \tag{26}$$

The study of the types of interaction potentials and the calculation of the approximate values of the force ranges, the kinetic energies, the internal energies, and the energies of activation make it possible to correct the minimal numbers Ni of

How to Use the Monte Carlo Simulation Technique? Application: A Study of the Gas Phase…

and the sizes (Lx, Ly, Lz) have to be discussed for statistical calculations.

For this MCS, the numerical chosen values are in Table 2.

; Qnp9 = 131).

As an example for our MCS calculation, we have:

The ratios Ns,i/Ns, H2 are too small (Table 3).

5.2.4 Calculation of statistical properties and some results of the calculations

between particles of the same species in the cell.

Let kp be the number of a species, kp = 1,…, 9. The minimal numbers Qnp(kp)

For numerical programming, according to the programming language used and according to the size (or the computational capacity) of the computer, it is necessary to find a judicious choice of the tables of integer or real values and which values would be useful to save all during simulation. Let Ncol,m be the maximum number of elastic collisions per particle, and let Ncycle be the number of cycles to average the

For radicals (e.g., SiH3), particle numbers Qnp(k) are very small; we take Qnp (k) = 10. These numbers cannot take value 1 or 0, even if a species k is in trace form in the gas. The value 0 for a species k means that any other species k' does not make a collision with the species k; and the value 1 means that we have no collisions

Qnp1, Qnp5, and Qnp9 are calculated from the volume of cell, the pressure, the temperature, and the total number of particles in the cell (Qnp1 = 0.81187824 \* 109

As we have chosen a stationary regime, we must reach the values and properties at equilibrium. The results of the simulation show this trend. In MCS, averaged values, distribution functions, autocorrelation functions, and correlation functions can be calculated. To ensure rapid convergence of calculations, it would be useful to look for statistically symmetric (or stationary or unsteady) parameters [26, 50].

• The number of Si2H6, SiH, and Si2H5 particles reaching the surface is negligible.

• Let Ns,i and Ns, H2 be the densities of a species i and H2 reaching the surface.

Lx (m) 4.68 10<sup>6</sup> 1 Qnp(SiH4) Qnp1 Ly (m) 4.68 10<sup>6</sup> 2 Qnp(SiH3) 10 Lz (m) 20.0 10<sup>3</sup> 3 Qnp(SiH2) 10

Ncol,m 500 5 Qnp(H2) Qnp5 Ncycle internal cycle 2000 6 Qnp(Si2H6) 10 Ncycle external cycle 200,000 7 Qnp(SiH) 10

Used quantities and parameters in calculations for the gas temperature Tg = 520 K.

Number of species Kp Initial number of particles in

4 Qnp(H) 10

8 Qnp(Si2H5) 10 9 Qnp(e) Qnp9

cell

;

particles and the sizes (Lx, Ly, Lz) of the elementary cell.

DOI: http://dx.doi.org/10.5772/intechopen.88559

simulation calculations.

Qnp5 = 0.20296956 \* 10<sup>9</sup>

Cell dimensions and steps for

collisions

Table 2.

131

where the pre-exponential factor is assumed to be the collision frequency factor and Kreac is the rate constant of the gas phase reaction.

The two colliding particles (e.g., the electron and SiH4 molecule) can interact by several reactions (R1, R2, R3, and R4 in Table 1); we choose randomly one of gas phase reactions occurring according to a, nonuniform discrete distribution reaction probability Prreac (i,j):

$$\text{Pr}\_{\text{reac}}(i, j) = \frac{K\_{\text{reac}}(i, j)}{\sum K\_{\text{reac}}(i, j)} \tag{27}$$

where <sup>P</sup>Kreacð Þ <sup>i</sup>; <sup>j</sup> is the sum of all rate constants of possible reactions between <sup>i</sup> and j.

All chemical systems go naturally toward states of minimum Gibbs free energy [21, 24]. A chemical reaction tends to occur in the direction of lower Gibbs free energy. To determine the direction of the reaction that is taking place, we use the old and new values of Kreac and the equilibrium constant with reactants and product concentrations. Each set of binary collisions can be related or converted into time. As cited in section (a), Table 1 gives gas phase reactions and corresponding rate constants used in this MCS.

To continue the simulation, after the elastic collision, particle i takes new values of components velocity and new mean free path; mean free path is taken from a normal (nonuniform) distribution (Gaussian distribution). If the collision is inelastic, we have to take a new particle.

From Metropolis algorithm, the scheme of this MCS is as follows:


#### 5.2.3 The choice of simulation parameters

Once the species are selected for the simulation model, an estimate of species densities should be made. Following the model of interaction and collisions between particles (binary, collective, etc.), a first choice of the minimum number Ni of particles of each species is made. A first estimate of the sizes (Lx, Ly, Lz) of the elementary cell is made.

How to Use the Monte Carlo Simulation Technique? Application: A Study of the Gas Phase… DOI: http://dx.doi.org/10.5772/intechopen.88559

The study of the types of interaction potentials and the calculation of the approximate values of the force ranges, the kinetic energies, the internal energies, and the energies of activation make it possible to correct the minimal numbers Ni of particles and the sizes (Lx, Ly, Lz) of the elementary cell.

Let kp be the number of a species, kp = 1,…, 9. The minimal numbers Qnp(kp) and the sizes (Lx, Ly, Lz) have to be discussed for statistical calculations.

For numerical programming, according to the programming language used and according to the size (or the computational capacity) of the computer, it is necessary to find a judicious choice of the tables of integer or real values and which values would be useful to save all during simulation. Let Ncol,m be the maximum number of elastic collisions per particle, and let Ncycle be the number of cycles to average the simulation calculations.

For this MCS, the numerical chosen values are in Table 2.

For radicals (e.g., SiH3), particle numbers Qnp(k) are very small; we take Qnp (k) = 10. These numbers cannot take value 1 or 0, even if a species k is in trace form in the gas. The value 0 for a species k means that any other species k' does not make a collision with the species k; and the value 1 means that we have no collisions between particles of the same species in the cell.

Qnp1, Qnp5, and Qnp9 are calculated from the volume of cell, the pressure, the temperature, and the total number of particles in the cell (Qnp1 = 0.81187824 \* 109 ; Qnp5 = 0.20296956 \* 10<sup>9</sup> ; Qnp9 = 131).

#### 5.2.4 Calculation of statistical properties and some results of the calculations

As we have chosen a stationary regime, we must reach the values and properties at equilibrium. The results of the simulation show this trend. In MCS, averaged values, distribution functions, autocorrelation functions, and correlation functions can be calculated. To ensure rapid convergence of calculations, it would be useful to look for statistically symmetric (or stationary or unsteady) parameters [26, 50].

As an example for our MCS calculation, we have:



Table 2. Used quantities and parameters in calculations for the gas temperature Tg = 520 K.

gas phase reaction. The difference between the kinetic energy of the two particles and

Theory, Application, and Implementation of Monte Carlo Method in Science and Technology

Ea ¼ �kBTln Kreac=νij

Prreacð Þ¼ <sup>i</sup>; <sup>j</sup> Kreacð Þ <sup>i</sup>; <sup>j</sup>

where <sup>P</sup>Kreacð Þ <sup>i</sup>; <sup>j</sup> is the sum of all rate constants of possible reactions between <sup>i</sup>

All chemical systems go naturally toward states of minimum Gibbs free energy [21, 24]. A chemical reaction tends to occur in the direction of lower Gibbs free energy. To determine the direction of the reaction that is taking place, we use the old and new values of Kreac and the equilibrium constant with reactants and product concentrations. Each set of binary collisions can be related or converted into time. As cited in section (a), Table 1 gives gas phase reactions and corresponding rate

To continue the simulation, after the elastic collision, particle i takes new values of components velocity and new mean free path; mean free path is taken from a normal (nonuniform) distribution (Gaussian distribution). If the collision is inelastic,

a. Choices of particle of spice i with random position, velocity, and mean free path; periodic boundary conditions are used to keep particles in the

c. Study of collision type (elastic, inelastic). If the collision is elastic the particle i move with a new velocity and mean free path, and we return to step (b). If the collision is inelastic particles i and j give new particles i' and j', according to Metropolis scheme, and we return to step (a) or (b). Periodic boundary conditions are used to keep particles in the elementary cell.

Once the species are selected for the simulation model, an estimate of species densities should be made. Following the model of interaction and collisions between particles (binary, collective, etc.), a first choice of the minimum number Ni of particles of each species is made. A first estimate of the sizes (Lx, Ly, Lz) of the

From Metropolis algorithm, the scheme of this MCS is as follows:

where the pre-exponential factor is assumed to be the collision frequency factor

The two colliding particles (e.g., the electron and SiH4 molecule) can interact by several reactions (R1, R2, R3, and R4 in Table 1); we choose randomly one of gas phase reactions occurring according to a, nonuniform discrete distribution reaction

� � (26)

<sup>P</sup>Kreacð Þ <sup>i</sup>; <sup>j</sup> (27)

the activation energy define the kind of collision (effective or not effective).

The activation energy is given by:

probability Prreac (i,j):

constants used in this MCS.

we have to take a new particle.

elementary cell.

b. Choices of random collision with a spice j.

d. At each step, we can note the different statistics.

5.2.3 The choice of simulation parameters

elementary cell is made.

130

and j.

and Kreac is the rate constant of the gas phase reaction.


collisions; it also makes it possible to make the statistics of the interactions with the

How to Use the Monte Carlo Simulation Technique? Application: A Study of the Gas Phase…

Other questions may be asked to account for molecular ions, surface and volume kinetics, or thin film formation. The techniques and different models of the MCS

The interconnection of the MCS with other models (MDS, hybrid model, fluid model, electromagnetic model, etc.) would allow answering more questions. The

surface. The results are compatible with [39, 51, 52].

DOI: http://dx.doi.org/10.5772/intechopen.88559

Author details

133

Fethi Khelfaoui\* and Oumelkheir Babahani

provided the original work is properly cited.

khelfaoui.fe@univ-ouargla.dz

Sciences, Kasdi Merbah Ouargla University, Algeria

\*Address all correspondence to: fethi.khelfaoui@gmail.com;

Physics Department, RPPS Laboratory, Faculty of Mathematics and Matter

© 2019 The Author(s). Licensee IntechOpen. This chapter is distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/ by/3.0), which permits unrestricted use, distribution, and reproduction in any medium,

(PP-MCS, MCS-PIC, kMCS) allow taking care of these questions.

methods can be applied to other specialties than the physical sciences.

Table 3.

Ratios Ns,i/Ns, H2 of particles reaching the surface compared to H2.


#### Table 4.

Ratios Ns,i/Nv,i of particles reaching the surface compared to volume.

