Keisuke HAYABUSA*
Masashi SHIMOYAMA*
Kenji AMAYA**
Tsutomu NAKADA*
*
Technologies, R&D Division
**
Tokyo Institute of Technology
Recently electro copper plating has been widely used for LSI interconnections. Because of the miniaturization of interconnections, the film thickness of wafer edges becomes extremely high. This phenomenon is particularly noticeable and called the terminal effect. In a plating process, film thickness distribution must as even as possible. Since extensive experience or numerous experiments are necessary to determine conditions optimal for plating, it has become common to include numerical analysis in building plating processes. The authors have proved the validity of a method developed by comparing numerical analysis to experimental data. This method applies vector-type boundary conditions and the region division method to issues subject to the base film resistance of wafers, multi-electrodes, shielding plate, and high resistors.
Keywords: Multiple Domain, Boundary Element Method, Numerical Analysis, Electrochemistry, Electroplating
Recently electro copper plating has been widely used for LSI interconnections. The technology provides advantages in cost, productivity, and performance in filling microtrenches and micropores. On the other hand, miniaturization of interconnections has required the seed and barrier layers to become thinner, causing the resistance to increase. During plating, the increase in resistance causes the phenomenon called “terminal effect”, where the film on the wafer edge becomes extremely thick, making film thickness distribution on the wafer less uniform 1). This phenomenon may increase the time required for the subsequent polishing process and cause problems, such as variations in electrical pad sizes within a chip, requiring the film thickness distribution to become as uniform as possible in the plating process.
Film thickness distribution can be leveled by changing the shape of electrodes, plating bathes, or shielding plates, or by adjusting electric current or plating time, etc. Since extensive experience or numerous experiments are required to determine optimum conditions for plating, it has become common to include numerical analysis in building plating processes. While Obara et al. 2), 3) have developed a numerical analysis program based on the FEM (Finite Element Method), and M. Purcar et al. 4) a program based on the boundary element method, the authors have developed a boundary element method program applicable to issues subject to a) the base film resistance of wafers, b) multi-electrodes, c) shielding plate, d) high resistors and e) a variety of symmetries, applying vectortype boundary conditions to a) and b), the region division method to c) and d), and the boundary element method algorithm for a variety of symmetries to e) 5)-13). By applying common approaches to various required functions in this way, an efficient program has been realized. This paper discusses the a) to e) approaches, compares the results with the experimental values, and verifies the validity of the developed program.
Figure 1 shows the configuration of the electro copper plating equipment 14). The wafer and anode are located facing each other in the rectangular solid plating bath filled with plating solution. Between the wafer and anode, the shielding plate and the paddle are installed to level the film thickness distribution and to stir the plating solution respectively. On the wafer surface serving as a cathode, the seed layer and barrier layer (hereinafter collectively referred to as the base film) are applied to conduct current for plating deposition and to prevent copper diffusion respectively. On the wafer periphery, current absorbing terminals (hereinafter referred to as terminals) are provided, serving as electrodes where the current running into the wafer base film finally flows out. Specified direct current is applied from the external power supply to the above configuration to deposit plating on the wafer surface.
Fig. 1 Electro plating bath
Figure 2 shows the configuration of the electro copper plating equipment, including its schematic diagram and circuit diagram. If contact resistance between the wafer and terminals is ignored, with power supply voltage φE , base film resistance Rf , anode polarization resistance Ra , cathode polarization resistance on the wafer surface RW , and plating solution resistance Rs , the value of current running around the center of the wafer Iin , and the value of current running around the wafer periphery Iout are calculated by the following equations.
When the base film is thin, Rf
cannot be ignored and Iin
becomes smaller than Iout
, resulting in increase in film thickness around the wafer periphery. In contrast, by placing the shielding plate and extending the current path to increase Rs
, the proportion of Rf
to the entire denominator of the equations (1) and (2) is reduced, making the film thickness distribution more uniform 1), 15).
When conducting quantitative examination using an actual plating bath, the plating surface needs to be treated as two-dimensional and the plating solution as three-dimensional. It must also be noted that a theoretical solution cannot be obtained since anode polarization resistance Ra
and cathode polarization resistance RW
are both nonlinear functions dependent on current density or electric potential. In this study, we introduced numerical analysis approaches, applying the finite element method for wafers and the boundary element method for plating solution to obtain the film thickness distribution 9)-11).
Fig. 2 Circuit diagram for plating bath
Figure 3 shows the governing equations and boundary conditions for the plating simulation used for this study. First, the following Laplace equation is used for the governing equation for the plating solution.
If the equation (3) is integrated twice and rearranged, the following boundary integral equation is obtained.
k represents electric conductivity of the solution, p the source point, and q the observation point. φ(q) and i (q) represent electric potential and current density respectively at the point q. φ* (q, p) and i* (q, p) represent fundamental solutions, and are calculated by the following equations when the three-dimensional field is taken into consideration.
As mentioned above, the conductive base film is applied to the wafer. When the film is thin, its resistance needs to be considered. Therefore, the wafer base film is modeled using the following Poisson equation.
kf represents electric conductivity of the base film, i w the density vector of current flowing into the wafer surface, and φf the vector of electric potential distribution inside the wafer base film. Since current value Im is specified for each electrode in a plating process, the following condition is added to this analysis.
n e represents the number of electrodes outside the wafer, and the value Iw of current running into the wafer is obtained by the following equation.
Distributions of electric potential and current density in the plating bath and wafer are obtained by applying each boundary condition to simultaneous equations formulated by discretizing the equations (4) to (9). At first, the following equation is obtained by applying a finite element to the equation (4) to discretize.
In this equation, H and G represent the matrices determined by geometric and material conditions in the analysis field, and φ and i represent boundary value vectors of electric potential φ and current density i in each element. The table shows the boundary conditions of general electroplating consisting of anode boundaries Γ A , wafer surface boundariesΓ w , and insulation boundaries (shielding plate, plating bath wall, etc.)Γ N . The anode and wafer boundaries with metal surfaces follow the polarization curve φ=−f (i) with current density i (=x) unknown. The insulation boundary has zero current density and unknown electric potential. Vectors of these unknown values are collectively represented by x =｛xA ，xw ，xN ｝.
Fig. 3 Governing equation and boundary condition
| Boundary Condition | Potential | Current Density |
| Anode Γ A | −fA （xA ） | xA |
| Wafer Γ W | −fW （xW ） | xW |
| Insulation Γ N | xN | iN （≡0） |
As mentioned previously, plating analysis needs to deal with the specified current values for electrodes (Equation (8)) and resistance of the wafer base film (Equation (7)). Since the analysis program becomes large if dealing with the values separately, this study applies the vector-type boundary conditions which can deal with these conditions in an integrated form.
The equations (11) and (12) show the relations between the unknown value xV of an element group and boundary values φV and iV .
This section applies the finite element method as a numerical analysis approach to equation (7) as a governing equation of the wafer base film, and describes the process of formulating the method as the vector-type boundary condition mentioned previously. Since the wafer base film and plating film thicknesses are significantly small in relation to the wafer size, they are handled as a two-dimensional problem and the film thickness distribution is used as a parameter. As shown in Figure 3 and Figure 4, assuming that the electric potential φf is zero at the terminals Γ D on the wafer edge, the current if is zero on the wafer edge Γ N excluding the terminals, and the density −iw of current flowing into the wafer surface, the following simultaneous equation for the entire wafer can be solved.
K represents the coefficient matrix determined by coordinates of each node,φn
f
the electric potential vector of the wafer node, and iw
the density vector of current flowing into each element on the wafer. Here, the mesh on the wafer is the part shared by the boundary element model of the solution and the finite element model of the wafer. While the boundary element model of this analysis uses a finite element with a finite value within the element, the finite element model uses a primary element with a value at each node. To share both element models, the average of node values within the element of the finite element model is handled as the element value.
Next, we explain the method of using the vectortype element conditions to express behavior inside the wafer base film. Figure 4 shows how the current −iw
flows from the solution into the wafer and flows out of the terminals through the base film. The base film resistance is large enough to generate electric potential distribution. For this reason, assuming the boundary conditions Γ
D
, and Γ
N
mentioned previously, conduct a finite element analysis in the case that the current with 1 A/m2 density flows into an element j to obtain electric potential distribution φj across all the elements on the wafer.
Then, application of the same unit current density to other elements ij＋1 is repeated until all the elements (the number of elements nf ) on the wafer are processed. The electric potential values thus obtained are arranged as column vectors and create the matrix represented by the following equation.
Electric potential distribution φf inside the wafer can be expressed by overlapping the density of current flowing into the elements. The distribution is calculated by the following equation using the obtained Af and the density vector iw of current flowing into the wafer.
The following simultaneous equation is formulated with the electric potential φw on the plating solution side by adding the electric potential gap on the wafer surface shown by the dashed line in Figure 3.
If this equation is expressed using the same vector f type boundary conditions as used for the equations (11) and (12), the following equations are obtained.
Fig. 4 Wafer model
Specifying current values is one of the important conditions in an electroplating process. This section explains the method of specifying current values for the electrodes including wafers and anodes using the vector-type boundary conditions mentioned previously. The current value Im to be specified for an electrode m is calculated by the following equation using the area Sj and the current density ij at the element j composing the electrode.
nm represents the number of elements on the electrode. If selecting a certain element (any can be selected) and numbering it as 1 (one), the equation (20) is transformed into the following using the current density i1 of the element.
In the case of specifying a current value, offset voltage φE by the power supply becomes unknown. Electric potential φA is calculated by the following equation including polarization resistance −fA (iA ) on the electrode surface.
From the equations (21) and (22), the unknown vector under the vector-type boundary conditions, as shown in the following equation, becomes the type where the element 1 has an electric potential component and the other elements have a current density component.
From the above, the vector-type boundary conditions are calculated by the following equations.
If the wafer surface is regarded as a current value specifying boundary, the program will become complicated. Here electrodes excluding the wafer are regarded as current values specifying boundaries. Accordingly, when the number of electrodes including the wafer is m＋1, that of the current value specifying boundaries is m.
When the equations (18), (19), (24), and (25) are substituted into the equation (10), the following equation is obtained.
Equation (26) is a nonlinear simultaneous equation in relation to the unknown vector x={xN , xA , xw } T and is analyzed using Newton-Raphson’s method. Since Jacobians for H and G are required to use this method, a Jacobian equation is used in the case that vector-type boundary conditions are included. First, equation (26) is turned into the following using the residual vector r.
The following equation is obtained by partially differentiating the residual vector r with an unknown vector x.
When noting the items related to the electrode (corresponding to the second line of the vectors) in equation (26), the electric potential vector and current density vector are calculated as follows.
Jacobian J is a matrix with a component in offdiagonal terms.
In plating processes, shielding plates or high resistors may be installed in the plating bath to level the film thickness distribution. Abe et al. have analyzed the case of installing a thin shielding plate in the plating bath by using a special weight function instead of a fundamental solution 5). Aoki et al. have analyzed the case of electric conductivity changing slowly in soil by deriving a fundamental solution in the case of linear change of electric conductivity 7). Since the analysis program becomes complicated if dealing with these cases separately, this paper discusses the use of the region division method to deal with the above issues 16), 17).
The region division method solves issues by applying a boundary integral equation to each region in the system consisting of several homogeneous regions to satisfy conditions compatible between regions. As shown in Figure 5, when each region with electric conductivity k1
, or k2
is divided into the regions Ω1 and Ω2, a discretized boundary integral equation is obtained for each region.
The superscript 1 or 2 corresponds to each region and the subscript I to the interface between the regions Ω1 and Ω2, and the boundary values on the boundary Γ1 in the region Ω1 and on the boundary Γ2 in the region Ω2 are represented by φ 1 , i 1 and φ 2 , i 2 respectively. If it is assumed that the electrical potential φ is continuous on the boundary Γ I and that the current iI running through the boundary Γ I is consistent, the following equations are obtained.
In equation (34) the signs differ between i I 1 and i I 2 to signify inflow and outflow of current. In consideration of these conditions, the following simultaneous equation is formulated from the boundary integral equations.
A solution can be obtained using the same method as equation (10) mentioned previously 8) .
Fig. 5 Multiple domain method
When a thin plate structure such as a shielding plate is included, the potential distribution obtained from the source points located on both sides of the thin plate may become almost uniform, the matrix becomes singular and can not be analyzed 5). As shown in Figure 5, however, singularity may be avoided by setting a virtual boundary Γ
I
to the opening of the shielding plate to divide the region and applying equations (31) and (32) to regions 1 and 2 respectively. Electric conductivity in these regions is assumed to be consistent (k
1＝k
2) and equations (33) and (34) are applied to the boundary conditions on the virtual boundary Γ
I
.
On the other boundaries on the shielding plate, insulation i is assumed to be zero (0). On the other hand, the region division method can apply to the analysis of the case in which a high resistor is inserted between the wafer and anode 1), 15). For example, if region 2 in Figure 5 has high resistivity, k
1 becomes higher than k
2, then the analysis should be conducted assuming that the interface between the high resistor and plating solution represented as Γ
I
is live.
We conducted the simulation using bump plating to prove effectiveness of the proposed three-dimensional simulation method of film thickness. Using the wafer holder with 8 terminals shown in Figure 6 and considering the symmetric profile of the plating bath, one fourth of the wafer shown in Figure 7 was analyzed. The model consists of the wafer, anode, and insulation sections, and the wafer section consists of region A where plating is deposited and region B where the current is shielded by the photoresist (hereinafter referred to as the mask). While current does not flow into the masked region B, current flowing into the unmasked region A flows to the terminals through the base film under the mask. The resistance in this flow was considered. The analysis was conducted using an electric conductivity 52.7 S/m in the plating solution, a current of 0.92 A for two and one half hours, and the polarization curve in Figure 8. This analysis evaluated the change in in-plane uniformity (=film thickness max.- min./average thickness) depending on whether region B is masked or not. Figure 9 (a) shows the film thickness distribution with the mask and Figure 9 (b) without the mask. In Figure 9 (a) where current runs in the wafer periphery, without the mask, the maximum film thickness is 0.124 mm measured in the vicinity of the terminals, the minimum thickness is 0.106 mm measured in the middle between the terminals, and the in-plane uniformity is 15.9 %. In Figure 9 (b), in the masked region, the film thickness distribution is zero (0), the maximum film thickness is 0.127 mm, the minimum thickness is 0.108 mm, and the in-plane uniformity is 17.2 %. When excluding region B in Figure 9 (a) not used in the product, the maximum film thickness is 0.116 mm and the in-plane uniformity in (a) is 8.6 %. From these results, the simulation has revealed that the in-plane uniformity can be substantially improved by depositing plating on the wafer periphery without even masking there.
Fig. 6 Wafer holder with 8 terminals
Fig. 7 3D simulation model
Fig. 8 Polarization curve
Fig. 9 Distribution of thickness on chips
We compared the film thickness distributions resulting from the numerical simulation with the experimental values, using bump plating. To optimize calculation involving the plating bath in Figure 6, we used cylindrical geometry inscribed inside the plating bath as a rectangular solid for approximation. Considering that the wafer holder has eight terminals, we used the cake-like 1/16 area of the wafer shown in Figure 10 as an analytical model. To this issue of symmetries, we applied the boundary element method algorithm for a variety of symmetries developed by the authors 13). The boundary conditions and number of element divisions were defined for each line segment in the cross section shown in Figure 10 (a), and when the number of circumferential element divisions was set as five, the numbers of elements and nodes were 320 and 380 respectively. Assuming that the wafer specifications and other conditions were the same as in the previous section, and that the average current density was 500 A/m2 , the actual plating thickness was measured using a resistance film thickness gauge. Figure 11 (a) and (b) show the results scaled by dividing the film thickness distributions obtained by the numerical simulation and the experiment respectively by the average thickness value. The comparison of Figure 11 (a) and (b) shows that the film consistently tends to be thick near the wafer center and terminals in both the simulation and experimental results, although the thickness distribution by the experiment shows a slight variation. In addition, the film thickness distribution from the wafer center toward periphery shown in Figure 12 reveals that the numerical simulation result mostly agrees with the experimental result, proving high accuracy of the simulation.
Fig. 10 Simulation model in case of bump plating
Fig. 11 Distribution of the Cu thickness
Fig. 12 The Cu thickness on the OA
In this study, we developed the boundary element method program that allows us to apply the region division method and vector-type boundary conditions and realized the numerical simulation with various functions for electro copper plating. We proved:
（1）
We can realize a numerical simulation program with many functions without expanding the program size by applying the vector-type boundary conditions to the wafer base film resistance and multiple electrodes and applying the region division method to the shielding plate and high resistor.
（2）
In the comparison between the film thickness distributions obtained by the numerical simulation and by the measurement using a resistance film thickness gauge, both tendencies were consistent. The result proves the high accuracy of the developed program, allowing us to verify its validity.
1) Kanda, H., Mine, J., Yamamoto, S., Kurashina, K., Hayabusa, K., Susaki, A., Nakada, T., Development of Cu damascene electroplating machines for interconnects on highly resistive substrates, EBARA ENGINEERING REVIEW, No.222 (2009), pp.10-18 (in Japanese).
2) Obara, K., Hatase, H. and Sakurai, T., Calculation of current distribution and surface finishing, The journal of the Surface Finishing Society of Japan, Vol.47, No.9 (1996), pp.751–757 (in Japanese).
3) Obara, K., Equation of the Film Thickness of Electroplating Using Current Distribution Analysis, The journal of the Surface Finishing Society of Japan, Vol.61, No.5 (2010), pp.366–370 (in Japanese).
4) M., Purcar, C., Munteanu and V., Topa, 3D electrode shape change simulation in electroplating, Rev. Roum. Sci. Techn. -Electrotechn. et Energ., Vol.58, No.3, pp.256–257 (2013).
5) Abe, K., Amaya, K., Aoki, S., Potential field BEM analysis for thin structure by special weighting function, Transactions of the Japan Society of Mechanical Engineers, Series A (2004), Vol.70, No.689, pp.31–35 (in Japanese).
6) Amaya, K., Naruse, N., Aoki, S., Miyasaka, M., Efficient boundary element method for large scale corrosion problem, Transactions of the Japan Society of Mechanical Engineers, Series A (1999), Vol.65, No.635, pp.1439–1497 (in Japanese).
7) Aoki, S., Takahashi, T., Amaya, K., Urago, M., Nishikawa, A., Boundary element analysis of galvanic corrosion of long structure buried in soil with gradient electrical conductivity, Journal of Society Material Science, Vol.46, No.2 (1997), pp.172–175 (in Japanese).
8) Aoki, S., Amaya, K., Miyasaka, M., Boundary element analysis on corrosion problems (1998), Shokabo, pp.35–36 (in Japanese).
9) Aoki, S., Amaya, K., Takazawa, H., Miyasaka, M., Boundary element simulation of electroplating on silicon wafer, The Journal of the Surface Finishing Society of Japan, Vol.51, No.4 (2000) pp.425–430 (in Japanese).
10) Hayabusa, K., Takemura, T., Amaya, K., Abe, K., Development of a plating analysis system for cyclic symmetry models, Proceedings of the Mechanical Engineering Congress, Vol.1, No.05-1 (2005), pp.55–56 (in Japanese).
11) Hayabusa, K., Amaya, K., Development of a plating simulation method for 3-Dimensional models, Proceedings of M&M 2006, No.06-4 (2006), pp.409–410 (in Japanese).
12) Hayabusa, K., Takemura, T., Amaya, K., Development of the optimization technology in planarizing process on semiconductor plating, Proceedings of the Mechanical Engineering Congress, Vol.6, No.07-1 (2007), pp.139–140 (in Japanese).
13) Hayabusa, K., Amaya, K., Boundary element algorithm and data structure for variety of symmetries, Transactions of the Japan Society of Mechanical Engineers, Series A, Vol.75, No.749 (2008), pp.21–26 (in Japanese).
14) Kuriyama, F., Minami, Y. and Kimura, M., Bump plating system, EBARA ENGINEERING REVIEW, No.207 (2005), pp.37–38 (in Japanese).
15) Hodai, M., Nakada, T., Kanda, H., Cu interconnect plating system, EBARA ENGINEERING REVIEW, No.207 (2005), pp.30–33 (in Japanese).
16) Yuuki, R., Kisu, H, Elastic analysis with boundary element method (1987), Baifukan, pp.165–170 (in Japanese).
17) Yuuki, R., Cao, G., Accurate BEM elaststatic analysis for very slender body, BEM, Vol.7 (1990), pp.53–56 (in Japanese).
Under the Scenes of our Lives High-pressure pump - Role and Application -
50% capacity boiler feed pump(BFP）playing an active role in a 1 000 MW thermal power plant
Large-capacity, Ultrahigh-efficiency, High-pressure Pumps for Seawater RO Desalination Delivered to Carlsbad Desalination Plant in the U.S.
Streamlines in crossover passage and velocity distributions at inlet of the second-stage impeller (Left:original,Right:optimized)
Discussion Meeting Symposium Ebara research system - Cooperation between research and business to create a new future -
Discussion Meeting (Mr. HIYAMA, Mr. SOBUKAWA, Mr. GOTO）
Under the Scenes of our Lives Standard Pumps - Essential Part of our Everyday Lives -
Examples of standard pumps
Inquiry about Ebara Engineering Review
Inquiry Form