Next Article in Journal
Polyphosphazenes: Multifunctional, Biodegradable Vehicles for Drug and Gene Delivery
Next Article in Special Issue
Energetic and Entropic Contributions to the Landau–de Gennes Potential for Gay–Berne Models of Liquid Crystals
Previous Article in Journal
Biodegradable Poly(butylene succinate) Composites Reinforced by Cotton Fiber with Silane Coupling Agent
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Multiscale Modeling of Chemical Vapor Deposition (CVD) Apparatus: Simulations and Approximations

Ernst-Moritz-Arndt-University of Greifswald, Felix-Haussdorffstr. 6, Greifswald D-17487, Germany
Polymers 2013, 5(1), 142-160; https://doi.org/10.3390/polym5010142
Submission received: 27 November 2012 / Revised: 16 January 2013 / Accepted: 22 January 2013 / Published: 5 February 2013
(This article belongs to the Special Issue Multiscale Simulations in Soft Matter)

Abstract

:
We are motivated to compute delicate chemical vapor deposition (CVD) processes. Such processes are used to deposit thin films of metallic or ceramic materials, such as SiC or a mixture of SiC and TiC. For practical simulations and for studying the characteristics in the deposition area, we have to deal with delicate multiscale models. We propose a multiscale model based on two different software packages. The large scales are simulated with computational fluid dynamics (CFD) software based on the transportreaction model (or macroscopic model), and the small scales are simulated with ordinary differential equations (ODE) software based on the reactive precursor gas model (or microscopic model). Our contribution is to upscale the correlation of the underlying microscale species to the macroscopic model and reformulate the fast reaction model. We obtain a computable model and apply a standard CFD software code without losing the information of the fast processes. For the multiscale model, we present numerical results of a real-life deposition process.
Classification: MSC:
35K25; 35K20; 74S10; 70G65

1. Introduction

In recent years, chemical vapor deposition (CVD) processes have received important applications to metal plates. Metallic or ceramic materials, such as SiC or a mixture of SiC and TiC, can be deposited in thin layers to substitute for expensive full metal plates. Our contributions are to apply such delicate multiscale models for simulating the CVD processes and reduce such models with respect to upscaling ideas to less complex and computable models (see [1]). We report the simulation results of a chemical vapor deposition (CVD) process. Such processes are applied to deposit thin films onto metallic or ceramic materials (see [2]). In the last few years, there has been much investigation of the optimization of such deposition processes. An example are thin films based on low temperature and low pressure processes with a mixture of standard applications to SiC and TiC (see [3]). We concentrate on deposing SiC films, which are important, but delicate to model and optimize with regard to obtaining a homogeneous deposition rate. Such homogeneous layers are important to achieve a stable nanolayer. We present a mixed model for the transport and kinetic processes of the CVD process with Tetramethylsilane as the precursor gas in a low temperature and low pressure plasma. We take into account the multiscale model of a large spatialtime-scale for the transport model and a small time-scale for the kinetic model of the CVD process. The plasma is modeled by an underlying quasi-equilibrium and neutral gas, which retards the precursor molecules in the kinetic model.
We use two software packages:
  • The macroscopic model (a transportreaction model with systems of coupled partial and ordinary differential equations) is simulated by UG/RNT (see [4]).
  • The microscopic model (a kinetic model with ordinary differential equations) is simulated by MATLAB (see [5]).
The present paper is organized as follows. In Section 1 and Section 2, we present the physical and mathematical model. Next, we simplify and reduce the original model to another model. In Section 3, we present the analytical and numerical methods that will be applied and the analysis of the coupled model equations. The numerical experiments are given in Section 4. In the conclusion, which is given in Section 5, we summarize our results.

2. Mathematical Model

In the following, the models are for the simulation of transport problems in the CVD apparatus. One can consider two scales:
  • Macro-scale of transport and reactions of the continuous species (scale of the apparatus);
  • Micro-scale of transport and reactions of the discrete particles (kinetic processes or scale of the atoms).
Here, we discuss the macro-scale and analytically embed the microscale of the reaction processes. We will discuss the following multiscale model:
  • Reactiondiffusion equations (see [6] (far-field problems));
  • Reaction equations that are embedded in the macroscale (see [7] (kinetic problems)).
We consider macroscopic problems based on small Knudsen numbers, Kn 0 . 01 1 . 0 . The Knudsen number (Kn) is the ratio of the mean free path λ to the typical domain size, L. As kinetic problems, we only consider the macroscopic chemical reaction between the clusters of species (see [7]). For a first overview of the apparatus, the full geometry (far-field) of the CVD apparatus is shown in Figure 1. A detailed graph with the dimensions of the apparatus is presented in Section 4.2.
Figure 1. Far field of the parallel chemical vapor deposition (CVD) apparatus.
Figure 1. Far field of the parallel chemical vapor deposition (CVD) apparatus.
Polymers 05 00142 g001
We consider the interesting deposition areas (near-field) in the apparatus, shown in Figure 2.
Figure 2. Near-field of the deposition area.
Figure 2. Near-field of the deposition area.
Polymers 05 00142 g002

2.1. Macroscopic Model for the Transport and Reaction Part

When gas transport is physically more complex due to combined flows in three dimensions, the fundamental equations of fluid dynamics become the starting point of the analysis. For our models with small Knudsen numbers, we can assume a continuum flow. The fluid equations can be treated with a NavierStokes or especially with a convectiondiffusion equation. Three basic equations, describing the conservation of mass, momentum and energy, are sufficient to describe the gas transport in the reactors (see [2]).
  • Continuitythe conservation of mass requires the net rate of mass accumulation in a region to be equal to the difference between the inflow and outflow rates.
  • NavierStokesmomentum conservation requires the net rate of momentum accumulation in a region to be equal to the difference between the in- and out-rate of the momentum, plus the sum of the forces acting on the system.
  • Energythe rate of accumulation of internal and kinetic energy in a region is equal to the net rate of internal and kinetic energy by convection, plus the net rate of heat flow by conduction, minus the rate of work done by the fluid.
We will concentrate on the conservation of mass and assume that energy and momentum are conserved (see [6,8]). Therefore, the continuum flow can be described as a convectiondiffusion equation:
( ϕ + ( 1 ϕ ) ρ R i ) t c i + · ( v c i D e ( i ) c i ) = λ i ( ϕ + ( 1 ϕ ) ρ R i ) c i + k = k ( i ) λ k ( ϕ + ( 1 ϕ ) ρ R k ) c k + Q ˜ i
where we have the following parameters:
ϕ : effective porosity (-); c i : concentration of the ith species, e.g., Si, Ti, C phase ( mol / mm 3 ) , v : velocity in the underlying plasma atmosphere ( mm / s ) , D e ( i ) : element specific diffusion-dispersion tensor ( mm 2 / s ) , λ i : decay constant of the ith species ( 1 / s ) , Q ˜ i : source term of the ith species [ mol / ( mm 3 s ) ] , k ( i ) : indices of the predecessors reactants of the ith species ( ) , R i : retardation factor due to plasma for theith species ( mm 3 / mol ) , ρ : Density of the plasma in the apparatus ( mol / mm 3 ) ,
where i = 1 , , M and M denote the number of species.
The effective porosity is denoted by ϕ and signifies the ratio of air to plasma in the apparatus environment. It says how much ionized plasma is filled with respect to the neutral gas (air). The transport term is indicated by the velocity v that presents the direction and the absolute value of the plasma flux in the apparatus. The velocity field is divergence-free. The kinetic constant of the ith species is denoted by λi. Hence, k(i) signifies the predecessor reactant species of the ith species, i.e., i consists of the results of the k-th species. The initial value is ci0, and we assume a Dirichlet boundary condition ci1(x, t) that is sufficiently smooth.

2.2. Microscopic Model for the Reaction Part

The kinetic processes involve reactions with different precursor gases. Such chemical reactions are derived in the work of Zhang and Huettinger (see [9]; for the available data, see [10]). In the following, we present the underlying reaction equations in the microscopic scales. These will later be embedded in equations at a macroscopic scale. The precursor of SiC is Tetramethylsilane, and we have the following reaction mechanism:
Si ( CH 3 ) 4 · CH 3 + · Si ( CH 3 ) 3
· CH 3 + H 2 CH 4 + · H
· CH 3 + Si ( CH 3 ) 4 CH 4 + ( CH 3 ) 3 Si C ˙ H 2
· Si ( CH 3 ) 3 + H 2 HSi ( CH 3 ) 3 + H ˙
2 · Si ( CH 3 ) 3 Si ( CH 3 ) 4 + Si ¨ ( CH 3 ) 2
Si ( CH 3 ) 4 + H 2 HSi ( CH 3 ) 3 + CH 4
HSi ( CH 3 ) 3 Si ¨ ( CH 3 ) 2 + CH 4
Si ¨ ( CH 3 ) 2 · Si ¨ CH 3 + · CH 3
The last reaction produces the deposition of the SiC.

2.3. Upscaling the Microscopic Model and Coupling to the Macroscopic Model

In the following, the upscaling of the microscopic model is done with respect to embedding the fast scales into the next coarser scales (see also Figure 3). The fast microscale reactions of the full reaction Equations (2)–(9) can be reduced with respect to embedding the intermediate reactions of (4)–(7). We then obtain the following reduced equation system:
Si ( CH 3 ) 4 · Si ( CH 3 ) 3 + · CH 3
· Si ( CH 3 ) 3 Si ¨ ( CH 3 ) 2 + · CH 3
Si ¨ ( CH 3 ) 2 · Si ¨ CH 3 + · CH 3
where we assume the transformation of Equation (8) with a fast reaction of H and given as:
HSi ( CH 3 ) 3 Si ¨ ( CH 3 ) 2 + CH 4 Si ¨ ( CH 3 ) 2 · Si ¨ CH 3 + · CH 3
Further, the reduced reaction Equations (10)–(12) can be upscaled by assuming a fast reaction of the · CH 3 species. We obtain the upscaled reaction equation:
2 · Si ¨ CH 3 SiC + CH 4 + Si + H 2
This Equation (14) can be applied to the macroscopic model (1). Therefore, we obtain an efficient computable model, where we can apply the macroscopic time steps.
The last reaction produces the deposition of the SiC.
Figure 3. Upscaling of the Microscopic model.
Figure 3. Upscaling of the Microscopic model.
Polymers 05 00142 g003

3. Analytical Methods and Numerical Methods

This section treats the underlying analytical methods and numerical methods to solve the multiscale models for the transportreaction Equation (2).

3.1. Multiscale Expansion (Embedding the Fast Scales)

We consider the multiscale equation
d y d t = Λ 1 y + 1 ϵ Λ 2 y
y ( 0 ) = y 0
where
y ( t ) = y 1 ( t ) y 1 ( t ) y m ( t ) , Λ k = λ 11 , k λ 12 , k λ 1 m , k λ 21 , k λ 22 , k λ 2 m , k 0 0 λ m 1 , k λ m 2 , k λ m m , k
where k = 1 , 2 , λ i j , k I R + with ϵ < < min i , j = 1 m λ i j , k and λ i j , k O ( 1 ) for all i , j { 1 , , m } , k = 1 , 2 .
Further, m is the number of reactants or species. If we omit the fast term 1 ϵ Λ 2 , there is an analytical solution y ( t ) = exp ( Λ 1 t ) y 0 , which is accurate on the slow scale with the O ( 1 ) time-scale, but has O ( 1 ) errors on the fast time-scale O ( ϵ ) .
To analyze the behavior on both scales, we have to extend the time scales:
y ( t ) y 0 ( t , τ ) + 1 ϵ y 1 ( t , τ ) + 1 ϵ 2 y 2 ( t , τ ) + 1 ϵ I y I ( t , τ )
where I I N + and τ = t ϵ is the slow time variable (see also [11]). We obtain equations ordered with respect to different scales. Physically, the multiscale solutions y 1 , y 2 , , are higher resolutions of the different fine scales (multiscales). While y 0 only resolves the scale O ( 1 ) , we obtain a resolution of finer scales with the higher orders: y 1 resolves O ( ϵ ) , y 2 resolves O ( ϵ 2 ) , and so on.
Substituting (18) into (15) yields:
y 0 t + Λ 1 y 0 + 1 ϵ y 1 t + Λ 1 y 1 + y 0 τ Λ 2 y 0 + 0
where d y i d t = y i t + ϵ y i τ | τ = 1 ϵ t with i = 1 , , I .
We thus obtain a system of equations that is valid to as high an order as possible when τ = 1 ϵ t . The leading order equation given by O ( 1 ) terms is:
y 0 t = Λ 1 y 0
The next order equation given by the O ( ϵ ) terms is:
y 1 t = Λ 1 y 1 + Λ 2 y 0 y 0 τ
For the fast scale-dependent derivative, we can assume for the quasi-stationary in the fast scale derivative:
y 0 τ 0
and we obtain:
y 1 t = Λ 1 y 1 + Λ 2 y 0
Such a multiscale equation can be solved analytically and can be embedded into the macroscopic transportreaction equations.

3.2. Analytical Solution of the Multiscale Reaction Equation

The multiscale reaction equation is an ordinary differential equation:
t R i c i = λ i R i c i + λ i 1 R i 1 c i 1
where i = 1 , , m , and we put λ 0 = 0 . The decay factors are λ i 0 . 0 , and the retardation factors are R i > 0 . 0 . The initial conditions are c 1 ( x , t 0 ) = c 01 and c i ( x , t 0 ) = 0 for i = 2 , , m .
We can solve these equations (cf. [12]):
c i = c 01 R 1 R i Λ i j = 1 i Λ j , i exp ( λ j t )
where i = 1 , , m . The solutions are defined for the case λ j λ k with j k and j , k 1 , , M . The factors Λ i and Λ j , i are given by:
Λ i = j = 1 i 1 λ j , Λ j , i = j = 1 j k i 1 λ k λ j
For equal reaction factors, we have derived the solution in ([13]). In the next subsection, we introduce the discretization of the diffusion-dispersion equation.

3.3. Numerical Methods for the TransportReaction Equation

For the numerical methods, we use finite volume methods for the space discretization (see [14,15]), and for the time discretization, we apply first order explicit or implicit Euler methods and second order CrankNikolson (CN) methods. For accurate results, we choose the second order CN method and accept the longer computational times, which are then needed. For fast and less accurate results, we can apply the cheaper explicit or implicit Euler methods. To embed the multiscale reaction equations, we use Godunov’s method for the multidimensional finite volume methods (see cf. [16]), and we could use one-dimensional analytical solutions of the convectionreaction equations.

3.4. Multiscale Embedding of the Reaction Parts into the Convection Part

To couple the upscaled microscopic reaction equation (23) with the macroscopic transport part, we apply Godunov’s method for discretization (cf. [16]). The formulation with the analytical solutions of the convection equations is extended to convectionreaction equations, while the multiscale expanded reaction equations can be used. We reduce the multi-dimensional equation to one-dimensional equations and solve each equation exactly. The one-dimensional solution is multiplied by the underlying volume, and we get the mass-formulation. The one-dimensional mass is embedded into the multi-dimensional mass formulation, and we obtain the discretization of the multi-dimensional equation.
The algorithm is as follows:
t c l + · v l c l = λ l c l + λ l 1 c l 1 with l = 1 , , m
The velocity vector v is divided by R l . The initial conditions are given by c 1 0 = c 1 ( x , 0 ) , or c l 0 = 0 for l = 2 , , m and the boundary conditions are trivial c l = 0 for l = 1 , , m .
We first calculate the maximal time-step for cell j and concentration i with the use of the total outflow fluxes:
τ i , j = V j R i ν j , ν j = k o u t ( j ) v j k
We get the restricted time-step with the local time-steps of cells and their components:
τ n min i = 1 , , m j = 1 , , I τ i , j
The velocity of the discrete equation is given by:
v i , j = 1 τ i , j
We calculate the analytical solution of the mass (cf. [13]) and we get:
m i , j k , o u t n = m i , o u t ( a , b , τ n , v 1 , j , , v i , j , R 1 , , R i , λ 1 , , λ i ) , m i , j , r e s t n = m i , j n f ( τ n , v 1 , j , , v i , j , R 1 , , R i , λ 1 , , λ i )
where a = V j R i ( c i , j k n c i , j k n ) , b = V j R i c i , j k n and m i , j n = V j R i c i , j n . Furthermore, c i , j k n is the concentration at the inflow boundary of cell j, and c i , j k n is the concentration at its outflow boundary.
The discretization with the embedded analytical mass is:
m i , j n + 1 m i , r e s t n = k o u t ( j ) v j k ν j m i , j k , o u t + l i n ( j ) v l j ν l m i , l j , o u t
where v j k ν j is the re-transformation for the total mass m i , j k , o u t of the partial mass m i , j k . In the next time-step, the mass is given by m i , j n + 1 = V j c i , j n + 1 , and in the old time-step, it is the rest mass for the concentration i. The proof is provided in ([13]).
In the next section, we derive an analytical solution for the benchmark problem (cf. [17,18]).

3.5. Discretization of the Source Terms

The source terms are part of the convectiondiffusion equations and are given as follows:
t c i ( x , t ) + v · c i D c i = q i ( x , t )
where i = 1 , , m , v is the velocity, D is the diffusion tensor, and the qi(x, t) are the source functions, which can be pointwise, linear in the domain.
The point sources are:
q i ( t ) = q s , i T t T , 0 t > T , , w i t h T q i ( t ) d t = q s , i
where q s , i is the concentration of species i at the source point x s o u r c e , i Ω over the whole time interval.
The line and area sources are:
q i ( x , t ) = q s , i T | Ω s o u r c e , i | , t T a n d x Ω s o u r c e , i 0 , t > T , with Ω s o u r c e , i T q i ( x , t ) d t d x = q s , i
where q s , i is the source concentration of species i at the line or area of the source over the whole time interval.
For the finite volume discretization, we have to compute:
Ω s o u r c e , i , j q i ( x , t ) d x = Γ s o u r c e , i , j n · ( v c i D c i ) d γ
where Γ s o u r c e , i , j is the boundary of the finite-volume cell Ω s o u r c e , i , j , which is a source area. We have j Ω s o u r c e , i , j = Ω s o u r c e , i , where j I s o u r c e , where I s o u r c e is the set of the finite-volume cells that includes the area of the source. The right-hand side of Equation (30) is also called the flux of the sources [19].
In the next subsection, we introduce the discretization of the diffusion-dispersion equation.

3.6. Discretization of the DiffusionDispersion Equation

We discretize the diffusiondispersion equation with implicit time discretization and the finite volume method for the equation:
t R c · ( D c ) = 0
where c = c ( x , t ) with x Ω and t 0 . The diffusiondispersion tensor D = D ( x , v ) is given by the Scheidegger approach ([20]). The velocity is v. The retardation factor is R > 0:0. The boundary values are denoted by n · Dc(x, t) = 0, where x ∈ Γ is the boundary Γ = ∂Ω, [21]. The initial conditions are given by c(x, 0) = c0(x).
We integrate (31) over space and time and obtain:
Ω j t n t n + 1 t R ( c ) d t d x = Ω j t n t n + 1 · ( D c ) d t d x
The time integration is done by the backwards Euler method, and the diffusion-dispersion term is lumped ([13]):
Ω j ( R ( c n + 1 ) R ( c n ) ) d x = τ n Ω j · ( D c n + 1 ) d x
Equation (33) is discretized over the space using Green’s formula.
Ω j ( R ( c n + 1 ) R ( c n ) ) d x = τ n Γ j D n · c n + 1 d γ
where Γ j is the boundary of the finite-volume cell Ω j . We use the approximation in space ([13]).
The spatial integration for Equation (34) is done using the mid-point rule over the finite boundaries and is:
V j R ( c j n + 1 ) V j R ( c j n ) = τ n e Λ j k Λ j e | Γ j k e | n j k e · D j k e c j k e , n + 1
where | Γ j k e | is the length of the boundary element Γ j k e . The gradients are calculated with the piecewise finite-element function ϕ l (see [22]), and we obtain:
c j k e , n + 1 = l Λ e c l n + 1 ϕ l ( x j k e )
With the difference notation, we get for the neighboring point j and l ([23]) and get the discretized equation:
V j R ( c j n + 1 ) V j R ( c j n ) = τ n e Λ j l Λ e { j } k Λ j e | Γ j k e | n j k e · D j k e ϕ l ( x j k e ) ( c j n + 1 c l n + 1 )
where j = 1 , , m .
In the next section, we discuss the numerical experiments.

4. Numerical Experiments

In the following, we present the numerical experiments of the microscale and macroscale simulations. An overview of the methods is given in Figure 4, where we present the different simulation methods applied for the microscale and macroscale simulations. The microscale simulation of the reaction equations gives the overview of the microscopic behavior with respect to the underlying temperature. The macroscale simulations of the transportreaction equations are compared with physical experiments. Here, we only applied upscaled simpler reaction parts, which embed the temperature characteristics in the macroscale equations. We could apply the physical results of the deposition rates and approximate our model equations with respect to the reaction and retardation parameters.
Figure 4. Simulation methods for the microscopic and macroscopic model.
Figure 4. Simulation methods for the microscopic and macroscopic model.
Polymers 05 00142 g004

4.1. Microscopic Experiment

We apply the full reaction equation (see also [10]) for the full kinetic equations of Tetramethylsilane (TMS) precursor (see Equations (2)–(9)). Based on the assumption of a fast reaction of the H species, we can apply the mesoscopic model equation, given as (see also Subsection 2.3):
Si ( CH 3 ) 4 · Si ( CH 3 ) 3 + · CH 3
· Si ( CH 3 ) 3 Si ¨ ( CH 3 ) 2 + · CH 3
Si ¨ ( CH 3 ) 2 · Si ¨ CH 3 + · CH 3
We assume the following velocity laws:
d [ S i ( C H 3 ) 4 ] d t = k × [ S i ( C H 3 ) 4 ]
d [ · S i ( C H 3 ) 3 ] d t = k × [ S i ( C H 3 ) 4 ] 0 . 9 k × [ · S i ( C H 3 ) 3 ]
d [ S i ¨ ( C H 3 ) 2 ] d t = 0 . 9 k × [ · S i ( C H 3 ) 3 ] 0 . 85 k × [ S i ¨ ( C H 3 ) 2 ]
d [ · S i ¨ ( C H 3 ) ] d t = 0 . 85 k × [ S i ¨ ( C H 3 ) 2 ]
where the temperature dependent reaction constant k is given by k ( T ) = 2 × 10 14 × exp [ 283 kJ mol 1 /( R T )] with R = 8.314472 J mol 1 K 1 . The derivation of the temperature dependent reaction constant k is discussed in the experimental work of [24,25]. The constants can be found in the NIST (National Institute of Standards and Technology) kinetics database ([26]).
In Figure 5, we show the differences between the different reaction temperatures, i.e., T = 573 K, 773 K, and 973 K, where we used the initial condition [Si(CH3)4]0 = 1 mol−1.
Figure 5. Decay of Si(CH3)4 (-), · S i ( C H 3 ) 3 ( · · ), : S i ( C H 3 ) 2 ( ), formation of · S i ¨ ( C H 3 ) ( · · ) and the summary of all concentrations (- -) at the temperatures 573 K (a), 773 K (b) and 973 K (c).
Figure 5. Decay of Si(CH3)4 (-), · S i ( C H 3 ) 3 ( · · ), : S i ( C H 3 ) 2 ( ), formation of · S i ¨ ( C H 3 ) ( · · ) and the summary of all concentrations (- -) at the temperatures 573 K (a), 773 K (b) and 973 K (c).
Polymers 05 00142 g005
Remark 1 The upscaled microscopic model shows the important influence of different temperatures. We obtain a slow reaction process at low temperatures and a fast reaction process at high temperatures. For the applied CVD process, an optimal temperature between 700 K and 900 K is appropriate. For such temperature regions, we see the dominance of the slow reaction rates: Such investigations allow application of our underlying multiscale reaction equation (24) for slow scales. Therefore, we can compute the macroscopic influence on the transport simulations, while we can upscale the microscopic scales to a simplified reaction process (see Subsection 2.3).

4.2. Test Experiment with SiC Deposition (Near-Field)

For all the experiments, we have the following parameters of the model, the discretization and the solver methods (Table 1).
Table 1. Physical and mathematical parameters.
Table 1. Physical and mathematical parameters.
Physical parameterMathematical parameter
Temperature, pressure, powervelocity, diffusion, reaction
T , p , WV , D , λ
In Figure 6, the underlying geometry of the apparatus is shown. The inflow of the precursor gases are at the left and right of the top of the apparatus, while the outflows are at the left and right bottom. The measured point ( 130 , 70 ) is in the middle of the deposition area at which the deposition rates could be measured.
Figure 6. The geometry of the apparatus with the measurement points (we apply ( mm ) as unit in the geometry).
Figure 6. The geometry of the apparatus with the measurement points (we apply ( mm ) as unit in the geometry).
Polymers 05 00142 g006

4.2.1. Parameters of the Model Equations

In the following, all the parameters of the model equation (2) are given in Table 2. Here, we have the physical experiments and approximate to the temperature parameters of T = 400, 600, 800 K. For the physical experiment, we have the following parameters (see Table 3).
Table 2. Model Parameters.
Table 2. Model Parameters.
density ρ = 0 . 5
mobile porosity ϕ = 0 . 333
diffusion D = 0 . 0
longitudinal dispersion α L = 0 . 0
transversal dispersion α T = 0 . 00
retardation factor R = 10 × 10 4 (Henry rate)
velocity field v = ( 0 . 0 , 4 . 0 × 10 8 ) t
decay rate of the species of 1st EX λ A B = 1 × 10 68
decay rate of the species of 2nd EX λ A B = 2 × 10 8 , λ B N N = 1 × 10 68
decay rate of the species of 3rd EX λ A B = 0 . 25 × 10 8 , λ C B = 0 . 5 × 10 8
Geometry (2d domain) Ω = [ 0 , 100 ] × [ 0 , 100 ] .
BoundaryNeumann boundary at
top, left and right boundaries.
outflow boundary
at the bottom boundary
Table 3. Approximated deposition rates and comparison to physical experiments.
Table 3. Approximated deposition rates and comparison to physical experiments.
WT P [ mbar ] R Si R C Physical ratio (Si:C)Numerical ratio (Si:C)
100700 9 . 7 × 10 2 4 × 10 4 2 × 10 4 0.5690.568
300700 9 . 7 × 10 2 2 . 3 × 10 4 2 × 10 4 0.7440.740
900700 9 . 7 × 10 2 1 . 35 × 10 4 2 × 10 4 0.9190.9
100400 1 × 10 1 2 × 10 4 0 . 7 × 10 4 0.6170.6103
500400 1 × 10 1 2 × 10 4 1 . 6 × 10 4 0.7570.745
500400 1 × 10 1 2 × 10 4 1 . 3 × 10 4 0.7040.691
900400 1 × 10 1 2 × 10 4 3 . 48 × 10 4 1.0101.017
900400 1 × 10 1 2 × 10 4 3 . 4 × 10 4 1.01.0
100400 4 . 5 × 10 2 4 . 7 × 10 4 0 . 1 × 10 4 0.3420.342
The discretization and solver method are the following:
  • For the spatial discretization method, we apply finite volume methods of the second order with the parameters in Table 4.
  • For the time discretization method, we used the CrankNicolson method (second order) with the parameters in Table 5.
  • The discretized equations are solved with the following methods, see the description in Table 6.
  • The initialization of sources of the equations are solved with the following parameters in Table 7.
Table 4. Spatial discretization parameters.
Table 4. Spatial discretization parameters.
spatial step size Δ x m i n = 1 . 56 , Δ x m a x = 2 . 21
refined levels6
limiterslope limiter
test functionslinear test function
reconstructed with neighbor gradients
Table 5. Time discretization parameters.
Table 5. Time discretization parameters.
Initial time-step Δ t i n i t = 5 10 7
controlled time-step Δ t m a x = 1 . 298 10 7 , Δ t m i n = 1 . 158 10 7
Number of time-steps 100 , 80 , 30 , 25
Time-step controltime steps are controlled with
the Courant number CFL m a x = 1
Table 6. Solver methods and their parameters.
Table 6. Solver methods and their parameters.
SolverBiCGstab (Bi-conjugate gradient method)
Preconditionergeometric multigrid method
SmootherGaussSeidel method as smoothers for
the multigrid method
Basic level0
Initial griduniform grid with 2 elements
Maximum Level6
Finest griduniform grid with 8192 elements
Table 7. Parameters of the source concentration.
Table 7. Parameters of the source concentration.
81 point sources of SiC at the position X = 10 , 11 , 12 , , 90 , Y = 20
Line source of H at the position x [ 5 , 95 ] , y [ 20 , 25 ]
Amount of the permanent source concentration ( · Si ¨ CH 3 ) source = 0.4, 0.7, 0.8, 0.85, 0.84, 0.82, 0.8,
0.6, 0.4, 0.2, 0.0., H source = 0 . 12
Number of time steps200

4.2.2. Numerical Results of the Model Equations

The numerical experiments now to be discussed are approximations to the SiC experiments. The underlying software tool is R 3 T , which was developed to solve discretized partial differential equations (see [4]). We use the tool to solve transportreaction equations. For the SiC, we obtain a different setup for the physical experiment, including the Bias voltage of the electric field, which is simulated as a retardation to the species. For the multiscale reaction equations, we can simplify the reaction process with respect to the slow scales. We consider an upscaled kinetic process, given by:
2 · Si ¨ CH 3 SiC + CH 4 + Si + H 2
In the following numerical experiment, we concentrate on the near-field computations of the deposition area (see Figure 2). We apply the transport-reaction parts (see Equation (1)) and the upscaled reaction (see Equation (45)).
We deal with the following parameters. Here we assume a constant velocity field and start with the species · Si ¨ CH 3 and H, which are given as point and line sources (see Table 8). We add some more H concentration to stabilize the scheme. We take here the concentration of · Si ¨ CH 3 as a point source, and the concentration of H is a line source. Further, we are interested in the relation between SiC and Si concentrations at the end of the deposition process. In Figure 7 and Figure 8, we present the concentration SiC, Si and H 2 , CH 4 after 100 and 200 time-steps. In the initialization, the amount of the Si C and Si species is not balanced; also, the amount of the H 2 species are too high. In such a situation, we would have a wrong deposition rate. In the later situation (see Figure 8), after 200 time steps, we see that the situation is balanced with respect to the SiC and Si concentrations. Here, fast reactions of · C H 3 and H have been passed, and we only have smooth transportreaction process. Now, the deposition of the layer is homogeneous and our rate is nearly 1 : 1 . In Figure 9, we show the results after the long deposition period of 200 time-steps. Here, the deposition rates are done with a 81 point sources experiment. Such a large amount of sources helps to homogenize the deposition in a large deposition region. We see a nearly constant deposition of the species SiC, while we dust small concentrations to the deposition area.
Table 8. Rate of the concentration.
Table 8. Rate of the concentration.
Rate at the end of the deposition at the deposited layer:
( · Si ¨ CH 3 ) s o u r c e , m a x : SiC t a r g e t , m a x
8 . 7 × 10 6 : 8 . 7 × 10 6 = 1 .
Figure 7. Experiment with moving point sources: SiC experiment after 100 time-steps, where a high concentration is red, a low concentration is blue (left figure: SiC concentration; middle figure: Si concentration; right figure: H concentration).
Figure 7. Experiment with moving point sources: SiC experiment after 100 time-steps, where a high concentration is red, a low concentration is blue (left figure: SiC concentration; middle figure: Si concentration; right figure: H concentration).
Polymers 05 00142 g007
Figure 8. Experiment with moving point sources: SiC experiment after 200 time-steps, where a high concentration is red, a low concentration is blue (left figure: SiC concentration; middle figure: Si concentration; right figure: H concentration).
Figure 8. Experiment with moving point sources: SiC experiment after 200 time-steps, where a high concentration is red, a low concentration is blue (left figure: SiC concentration; middle figure: Si concentration; right figure: H concentration).
Polymers 05 00142 g008
Figure 9. Deposition rates for the 81 point sources experiment (x-axis: time in 10 9 s, y-axis: concentration in mol / mm 3 ).
Figure 9. Deposition rates for the 81 point sources experiment (x-axis: time in 10 9 s, y-axis: concentration in mol / mm 3 ).
Polymers 05 00142 g009
Remark 2 The numerical experiments in the near-field can be approximations of the real-life physical experiments. Both experiments show the influence of temperature, while for low temperatures, we can assume we are dealing with slow time-scale reaction equations. In such regimes, we obtain the best results with multiple sources and long-time depositions. We apply further different experimental situations, and the best deposition result is obtained with low temperature and high power assumptions. At least homogeneous concentrations below the deposition area can be achieved with a large amount of sources. The near-field simulations obtain an optimum at the low temperature of 400 °C and a high plasma power of about 900 W. Such results are also obtained in our physical studies (see [27]).

5. Conclusions

We have presented a multiscale model for chemical vapor deposition processes. While for higher temperature regions, fast reaction rates are important, we embed such characteristics with multiscale expansions in our underlying transportreaction equations. In the real-life experiments, we see that only the slow reaction rates are important, because of the necessary low temperature regime to obtain an optimal homogeneous deposition. Approximations for the real-life experiments are made for a realistic apparatus with transport reactions.
The embedding of the multiscale reaction equations allows discretizing with a fast finite volume method and applying our underlying software code to the complex material functions of the model. We present models for the stoichiometry for SiC depositions and present their experiments. In the future, near- and far-field simulations will be able to derive the optimal parameter settings and be able to forecast the results of real-life experiments. Such simulations will then help to reduce the number of physical experiments that need to be carried out and give direction to future expensive physical experiments. In our future work, we will concentrate on further implementations of multiscale methods to higher temperature regimes.

References

  1. Turner, R. Computable Models, 2nd ed.; Springer: New York, NY, USA, 2009. [Google Scholar]
  2. Ohring, M. Materials Science of Thin Films, 2nd ed.; Academic Press: New York, NY, USA, 2002. [Google Scholar]
  3. Barsoum, M.W.; El-Raghy, T. Synthesis and characterization of a remarkable ceramic: Ti3SiC2. J. Am. Ceram. Soc. 1996, 79, 1953–1956. [Google Scholar] [CrossRef]
  4. Fein, E. Software Package r3t—Model for Transport and Retention in Porous Media; Technical Report GRS-192; Gesellschaft für Anlagen- und Reaktorsicherheit (GRS) mbH: Braunschweig, Germany, 2004. [Google Scholar]
  5. Downey, A.B. Physical Modeling in Matlab(R); CreateSpace Independent Publishing Platform: Los Angeles, CA, USA, 2009. [Google Scholar]
  6. Gobbert, M.; Ringhofer, C. An asymptotic analysis for a model of chemical vapor deposition on a microstructured surface. SIAM J. Appl. Math. 1998, 58, 737–752. [Google Scholar]
  7. Geiser, J.; Arab, M. Modelling for chemical vapor deposition: Meso- and microscale model. Int. J. Appl. Math. Mech. 2008, 4, 24–45. [Google Scholar]
  8. Geiser, J. Numerical simulation of a model for transport and reaction of radionuclides. In Proceedings of the Large Scale Scientific Computations of Engineering and Environmental Problems, Sozopol, Bulgaria, 6–10 June 2001.
  9. Zhang, W. CVD of SiC from methyltrichlorosilane. Part II: Composition of the gas phase and the deposit. Chem. Vap. Depos. 2001, 7, 173–181. [Google Scholar]
  10. Geiser, J.; Roehle, R. Kinetic processes and phase-transition of CVD processes for Ti2SiC3. J. Converg. Inf. Technol. 2010, 5, 9–32. [Google Scholar]
  11. Weinan, E. Principles of Multiscale Modelling; Cambridge University Press: Cambridge, UK, 2011. [Google Scholar]
  12. Bateman, H. The solution of a system of differential equations occurring in the theory of radioactive transformations. Proc. Camb. Philos. Soc. 1910, 15, 423–427. [Google Scholar]
  13. Geiser, J. Gekoppelte Diskretisierungsverfahren für Systeme von Konvektions-Dispersions-Diffusions-Reaktionsgleichungen. Ph.D. Thesis, Universität Heidelberg, Heidelberg, Germany, 2003. [Google Scholar]
  14. Geiser, J. Discretisation methods with analytical solutions for convection-diffusion-dispersion-reaction equations and applications. J. Eng. Math. 2006, 79, 1953–1956. [Google Scholar]
  15. Khattri, S.K. Analysing an adaptive finite volume for flow in highly heterogeneous porous medium. Int. J. Numer. Methods Heat Fluid Flow 2008, 18, 237–257. [Google Scholar]
  16. Leveque, R. Finite Volume Methods for Hyperbolic Problems; Cambridge University Press: Cambridge, UK, 2002. [Google Scholar]
  17. Higashi, K.; Pigford, T. Analytical models for migration of radionuclides in geologic sorbing media. J. Nucl. Sci. Technol. 1980, 17, 700–709. [Google Scholar] [CrossRef]
  18. Jury, W.; Roth, K. Transfer Functions and Solute Movement through Soil; Bikhäuser: Basel, Switzerland, 1990. [Google Scholar]
  19. Frolkovič, P. Flux-based Methods of Characteristics for Transport Problems in Groundwater Flows Induced by Sources and Sinks. In Computational Methods in Water Resources Volume II; Hassanizadeh, S.M., Schotting, R.J., Gray, W.G., Pinder, G.F., Eds.; Elsevier: Amsterdam, The Netherlands, 2002; pp. 979–986. [Google Scholar]
  20. Scheidegger, A. General theory of dispersion in porous media. J. Geophys. Res. 1961, 66, 32–73. [Google Scholar] [CrossRef]
  21. Frolkovič, P. Flux-based method of characteristics for contaminant transport in flowing groundwater. Comput. Vis. Sci. 2002, 5, 73–83. [Google Scholar]
  22. Cai, Z. On the finite voulme method. Numer. Math. 1991, 58, 7713–7735. [Google Scholar]
  23. Frolkovič, P.; de Schepper, H. Numerical modelling of convection dominated transport coupled with density driven flow in porous media. Adv. Water Resour. 2001, 24, 63–72. [Google Scholar] [CrossRef]
  24. Loumagne, F.; Langlais, F.; Naslain, R. Kinetic laws of the chemical process in the CVD of SiC ceramics from CH3SiC13–H2 precursor. J. Phys. IV 1993, 3, 527–533. [Google Scholar]
  25. Langlais, F.; Loumagne, F.; Lespiaux, D.; Schamm, S.; Naslain, R. Kinetic Processes in the CVD of SiC from CH3SiC13–H2 in a Vertical Hot-Wall Reactor. J. Phys. IV 1995, 5, 105–115. [Google Scholar]
  26. NIST (National Institute of Standards and Technology). Kinetics Database, Reaction: (CH3)4Si→·CH3+(CH3)3Si·. Available online: http://kinetics.nist.gov/kinetics/Detail?id=1972CLI/GOW53-61:1 (accessed on 28 January 2013).
  27. Geiser, J.; Arab, M. Model of PE-CVD apparatus: Verification and simulations. Math. Problems Eng. 2010. [Google Scholar] [CrossRef]

Share and Cite

MDPI and ACS Style

Geiser, J. Multiscale Modeling of Chemical Vapor Deposition (CVD) Apparatus: Simulations and Approximations. Polymers 2013, 5, 142-160. https://doi.org/10.3390/polym5010142

AMA Style

Geiser J. Multiscale Modeling of Chemical Vapor Deposition (CVD) Apparatus: Simulations and Approximations. Polymers. 2013; 5(1):142-160. https://doi.org/10.3390/polym5010142

Chicago/Turabian Style

Geiser, Juergen. 2013. "Multiscale Modeling of Chemical Vapor Deposition (CVD) Apparatus: Simulations and Approximations" Polymers 5, no. 1: 142-160. https://doi.org/10.3390/polym5010142

APA Style

Geiser, J. (2013). Multiscale Modeling of Chemical Vapor Deposition (CVD) Apparatus: Simulations and Approximations. Polymers, 5(1), 142-160. https://doi.org/10.3390/polym5010142

Article Metrics

Back to TopTop