The Implementation of the Surface Charging Effects in Three-Dimensional Simulations of SiO2 Etching Profile Evolution ()
Keywords:Plasma Etching; Level Set Method; Profile Charging; Finite Elements Method
1. Introduction
The profile surface evolution itself in plasma etching, deposition and lithography development is a significant challenge for numerical methods for interface tracking itself. Level set methods for evolving interfaces [1,2] are specially designed for profiles which can develop sharp corners, change topology and undergo orders of magnitude changes in speed. Additionally, etching submicron features in insulating materials such as SiO2 or materials masked with insulators, e.g., photoresist-masked polysilicon, are accompanied by feature charging effects caused by the differences in electron and ion angular distributions [3,4]. In general, the motion of ions incident at the wafer surface on the electrode is more directional than that of electrons. The difference causes the electron shading effect at the etched feature, thus generating positive potential at the deep trench bottom. This charging effect induces many serious plasma process induced damage problems such as bowing, trenching, reactive ion etching lag, and notching [5]. The reduction in the device size and multilayer structures requires a high-aspect ratio in the SiO2 etching. As the aspect ratio increases, the charging effect becomes more profound. The negative charge that might develop on the sidewalls can bend the trajectories of the ions passing close to them. This bending effect may contribute to micro trenching. In addition, etching rate reduction with time and notching of polysilicon sidewalls during the over etching step have been attributed in part to feature charging effects.
In this paper we present some results of our efforts to include charging effects in 3D level set etching profile simulator developed earlier. In the next section, we describe shortly the most important ingredients of the etching profile simulation procedure: the profile surface movement via level set method, particle fluxes calculations using Monte-Carlo method and the electric field (generated by the profile charge) calculations on the basis of the finite elements method.
2. Simulation Procedure
2.1. Level Set Method
The basic idea behind the level set method is to represent the surface in question at a certain time t as the zero level set (with respect to the space variables) of a certain function j (t, x), the so called level set function shown in Figure 1. The initial surface is given by {x|j (0, x) = 0}. The evolution of the surface in time is caused by the surface processes in the case of the etching. The velocity of the point on the surface normal to the surface will be denoted by R(t, x), and is called velocity function. The velocity function generally depends on the time and space variables and we assume that it is defined on the whole simulation domain. At a later time t > 0, the surface is as well the zero level set of the function j (t, x). Namely, it can be defined as a set of points {x|j (t, x) = 0}. This leads to the level set equation in Hamilton-Jacobi form [6]:
(1)
in the unknown function j (t, x), where Hamiltonian function is given by H = R (t, x)|Ñφ(t, x)| and where j (0, x) = 0 determines the initial surface. Several approaches for solving level set equations exist which increase accuracy while decreasing computational effort. The most important are narrow band level set method, widely used in etching process modeling tools, and recently developed sparse-filed method [3]. The sparse-field method use an approximation to the level set function that makes it feasible to recompute the neighborhood of the zero level set at each time step. As a result, the number of computations increases with the size of the surface, rather than with the resolution of the grid.
Figure 1. Examples of various level set function depending on the initial surface.
We say that Hamiltonian function H is convex if the following condition is fulfilled:
(2)
where 
is a partial derivative of φ(t,x) with respect of xi. The non-convex Hamiltonians are characteristic for plasma etching and deposition simulations. During these processes the etching (deposition) rate, that defines the surface velocity function R(t,x), depends on the geometric characteristics of the profile surface itself, or more precisely, on the angle of the incidence of the incoming particles. In the cases under study here we shall consider an etching beam coming down in the vertical direction. These conditions are characteristic for ion milling technology, but angular dependence of the etching rates appears, more or less, in all etching processes. The simplest finite difference scheme that can be applied in these cases is the Lax-Friedrichs scheme, one that relies on the central difference approximation to the numerical flux function, and preserves monotonicity through a secondorder linear smoothing term. It is shown [6] that it is possible to use the Lax-Friedrichs scheme in conjunction with the sparse field method, and to preserve sharp interfaces and corners by optimizing the amount of smoothing in it. This is of special importance in the simulations of the etching processes in which spatially localized effects appear, like notching and microtrenching.
2.2. Particle Fluxes Calculations
Our simulation concept is similar in spirit to the 2D simulations presented in [7,8] and especially [9], where charging effects in 3D rectangular trench were analyzed using finite difference method. Monte Carlo technique is the only feasible method for calculating particle fluxes in 3D geometries. Trench wall charging strongly influences the charged particles motion and, consequently, particle fluxes which themselves determine the local etching rates. The etching process in medium/high density fluorocarbon plasmas is believed to consist of concurrent etching (of the SiO2 substrate in our case) and deposition (of a fluorocarbon polymer layer) phenomena [10]. Here the deposition process is neglected for the sake of simplicity. So, we consider only a simple case of chemical etching of SiO2 under Ar+ ions bombardment. It is also assumed that the electrons are absorbed at the hitting points, while the neutrals can be absorbed or diffusively reflected, once or many times depending on its sticking coefficient. The positive ions can be absorbed, or secularly reflected, depending on their energy and incident angle. It is assumed that charged particles pass on their charge when they hit the surface, and that this charge does not move after that, what is reasonable for insulating materials. At the boundaries above the profile surface, periodic conditions are assumed.
The energy and angular distribution functions for Ar+ ions and electrons, used as the boundary conditions for the calculations of electron and ion fluxes incident on the profile surface, are obtained by particle-in-cell (PIC) calculations using XPDC1 code [11-13] and shown in Figure 2.
(a)
(b)
Figure 2. a) Electron Energy Distribution Function (EEDF) and b) Ion Energy Distribution Function (IEDF) calculated by using PIC/MCC code [11,12].
2.3. Electric Fields Calculations
Finite difference method cannot be used when the geometry of the problem is irregular, as it is usual for the etching profiles. So, it naturally to use finite elements method. The calculations of the electric fields in this paper are performed by integrating a general finite element solver GetDP [14,15] in our simulation framework. GetDP is a thorough implementation of discrete differential forms calculus, and uses mixed finite elements to discretize de Rham-type complexes in one, two and three dimensions. Meshing of the computational domain is carried out by TetGen tetrahedral mesh generator [16]. TetGen generates the boundary constrained high quality (Delaunay) meshes, suitable for numerical simulation using finite element and finite volume methods. As a part of the post-processing procedure the electric fields, obtained on the unstructured meshes, are recalculated on the Cartesian rectangular domains containing the regions of the particles movement. In this manner, the electric field on the particles could be calculated by simple trilinear interpolation. Details about the integration of the whole simulation cycle and the etching rate calculations using particle fluxes could be found elsewhere.
3. Results and Discussions
Here we will present some calculations illustrating our approach in etching profile evolution simulation. Although a complete self-consistent cycle that includes profile charging and its influence on the charged particles motion is already implemented, we do not have representative results yet. The reason is the lack of appropriate computational resources, since the Monte Carlo step requires enormous number of particles in every level set step, in order to get satisfying statistics. So, here we present some results of the profile charging calculations for the fixed profile surface.
All the calculations are performed on 128 × 128 × 384 rectangular grid. The initial profile surface is a rectangle deep with dimensions 0.1 × 0.15 × 0.1 mm. Above the profile surface is the trench region. From its top the particles involved in etching process come from, while bellow is the non-etched material. In the case of the ion enhanced chemical etching the dependence of the surface velocity on the incident angle is simple: V = V0 cosθ. In Figure 3 the evolution of the etching profile under this assumption is presented. This is the simplest form of angular dependence, but it describes the ion enhanced chemical etching process correctly [11]. In this case we expect that the horizontal surfaces move downward, while the vertical ones stay still. This figure shows that it with optimal amount of smoothing gives minimal rounding of sharp corners, while preserving the numerical stability of the calculations. In Figure 4 the surface mesh of
Figure 3. Etching profile evolution during equidistant time steps (0 s, 30 s, 60 s, 90 s and 120 s).
the tetrahedral mesh corresponding to the last profile in Figure 3 is shown. This mesh defines the profile surface for the subsequent particle flux and electric field calculations.
Collisions in the simulation domain are neglected because mean-free path of ions and electrons is much longer than the simulation domain size. Also, the probability that more than two particles exist in the simulation space at the same time is very low because the life time of ions or electrons in the simulation region is much shorter than the incident time interval. The particle fluxes on profile are given in Figure 5. The neutral and the electron fluxes distributions are similar: they are concentrated on the top of the wafer and upper portion of the sidewalls, while the ions, due its directionality, are accumulated at the bottom of the profile. Solving the Laplace equation with such a charge distribution as boundary condition leads to the potential cut depicted in Figure 6.
Figure 5. Space charge density on the profile surface [Q/m2].
Figure 6. Potential cut obtained by solving Laplace equation.
The maximal value of the potential at the bottom of the profile is somewhat bigger than the average ion energy, while the top surface is slightly negatively biased.
The electric field vectors on the profile surface are shown in Figure 7. It is obvious the field stimulate notch creation by directing the ions toward sidewalls. At the very the top the field accelerate the ions downward, but it quickly changes direction at the lower positions of the profile. The zero field border height is determined by the etching plasma parameters.
4. Conclusion
Plasma processing is indeed one of the crucial modern technologies that have made a significant impact on fabrication of integrated circuits (IC) and consequently in electronics and consumer goods in general. Damage to ICs during manufacturing as a result of charging of the dielectrics during finalization of interconnects is both reducing the profitability and reducing the ability to reach large sizes of microchips and make complex system integration on a single chip. So, realistic three dimensional etching profile simulations are urgently needed,
Figure 7. Electric field on the profile surface [V/m].
but still lacking. It is understandable that realistic calculations should include complex surface reactions set, the effects of the profile charging, polymer deposition, as well as better statistic (greater number of particles) in the Monte Carlo step of the calculations, which is usually limited by the available computational resources.
Acknowledgements
The present work has been supported by the Ministry of Education and Science, Republic of Serbai, Projects O171037 and III41011.