Subsections

• 4.4.1 Korhonen's Model
• 4.4.2 Analytical Solution of Korhonen's Equation
• 4.4.3 The Sink/Source Term $\gamma$
• 4.4.4 An Explanation of Black's Law
• 4.4.5 Sarychev's Model
• 4.4.6 Finite Element Discretization of the Basic Equation

4.4 Prediction of the Voids Nucleation Sites

In the case of copper interconnects there are several types of vacancy diffusivities of importance, namely, diffusion in the bulk material, grain boundary diffusion, the diffusion along the copper/barrier interfaces, and the diffusion along copper/cap-layer interface. The bulk diffusion coefficient can be neglected for temperatures well below the melting point [55]. The simpliest way to take into account all diffusivities is by introducing an effective diffusion coefficient dependent on the diffusivity coefficients of all paths and on microstructural parameters such as grain size and grain boundary thickness [77]. Since most copper grains in a Damascene lines are not larger than 0.3 $\mu$m the majority of the interconnects used today will inevitably have a bamboo structure, so that averaging in the sense of an effective diffusion coefficient is feasible. However, fast diffusivity paths along copper/barrier and copper/cap-layer interfaces will influence the vacancy dynamics strongly depending on the specific interconnect layout. This influence can not be incorporated with simple averaging terms. The diffusivity must be a cumulative value as used in [55]

$$D_v = D_b+D_{gb}\Bigl(\frac{\delta_{gb}}{d}\Bigr)+D_{Cu/b}q_{Cu/b}+D_{Cu/N}q_{Cu/N}.$$ (4.12)

$D_b$, $D_{gb}$, $D_{Cu/b}$, and $D_{Cu/N}$, represent the diffusivity through the bulk, along the grain boundaries, copper/barrier interfaces, and copper/caplayer interfaces, respectively. $\delta_{gb}$ is the width of the grain boundary and $d$ the average length of a grain boundary. Coefficients $q_{Cu/b}$ and $q_{Cu/N}$ depend only on the layout geometry.

To predict the location of the void nucleation sites and the time needed for a void to nucleate a multitude of microstructural, chemical, process-related, and environment-related variables have to be taken into account. The microstructural characteristics which need to be considered are the grain boundary structure and the orientation of the crystal axes inside the grains. Numerous experiments [55,8] have shown that the design of the interconnect, passivation, and barrier layers, together with the material composition of these layers, has a significant impact on the localization of void nucleation sites.
Of crucial importance for the void-nucleating condition are the dynamics of the crystal vacancies of the interconnect metal. Depending on crystal texture of the used copper interconnect and it's geometry there are several possible paths for the vacancy diffusion.

4.4.1 Korhonen's Model

The coupling of the vacancy dynamics with the stress development was first introduced by Korhonen in [51]. The analysis was carried out for a passivated straight aluminum interconnect line with a columnar grain structure with the grain boundaries perpendicular to the substrate. The author assumes that vacancy diffusion along grain boundaries equilibrates the stress in such a way that it reduces to a hydrostatic state on a time scale that is short compared to that for long range diffusion along the entire line length. The recombination and generation of vacancies changes the concentration of the available lattice sites, which consequently influences the hydrostatic stress distribution. Specifically, the loss of the available lattice sites yields an increase of the hydrostatic stress.

This assumption greatly simplifies the analysis und allows an analytical solution of the governing equation. However, this approach is eglible only for passivated straight interconnects with a specific aluminum crystal texture.

Metal self-diffusion occurs via a vacancy exchange mechanism. The atomic flux $\mathbf{J}_A$ is equal and opposite the flux of vacancies $\mathbf{J}_V$. Therefore, the electromigration induced vacancy flux resulting from an applied electric field $\mathbf{E}$ is given by [51],

$$\mathbf{J}_V=C_V\frac{Z^*eD_V}{kT}\mathbf{E},$$ (4.13)

and correspondingly, the atomic flux,

$$\mathbf{J}_A=- C_A\frac{Z^*eD_V}{kT}\mathbf{E},$$ (4.14)

where $C_A$, $D_V$ and $C_V$, $D_V$ are the concentration and diffusion coefficient of atoms and vacancies, respectively. $Z^*$ is the effective valence (Section 4.2).
In [51] is assumed that material transport along an interconnect line is carried out by grain boundary diffusion alone. This approach is justified in the case of an aluminum-based metalization when the width of the interconnect line is much larger than the average grain size. Then for the given boundary the vacancy diffusion coefficients $D^{GB}_V$ and the effective diffusion coefficient are related,

$$D_V\thickapprox \frac{ \delta D^{GB}_V}{d},$$ (4.15)

where $\delta$ is the grain boundary width and $d$ is the grain size. As it was shown by Herring [78,79], diffusive fluxes of atoms (vacancies) arise due to potential differences between different locations of the interconnect line and the chemical potential given by,

$$\mu = \mu_0 - \Omega \sigma,$$ (4.16)

where $\Omega$ is the atomic volume, and $\sigma$ is the tensile stress normal to the grain boundary or the interface.

One of the basic assumptions made by Korhonen [51], is that electromigration will deposit atoms in grain boundaries in such a way that mechanical stress very fast becames uniform (compared with the self-diffusion along the interconnect) at a any particular cross section of straight interconnect (perpendicular to the electrical current direction). This assumption significantly simplifies analysis because the mechanical stress gradient is non-zero only in the electrical current direction where it opposites electromigration.

The model transport equations can be expressed in terms of either vacancies or atoms but in most cases [51] vacancy are chosen. The net vacancy flux along the length of the interconnect line induced by the gradient of the chemical potential and electromigration, taking into account (4.16) can be written as,

$$J_V=-C_V\frac{D_V}{kT}\Bigl(\Omega \frac{\partial \sigma}{\partial x} - Z^*e E\Bigr).$$ (4.17)

The equilibrium vacancy concentration in the presence of mechanical stress is given by,

$$C_V=C_V^{eq} (\Omega\sigma/kT),$$ (4.18)

and the vacancy flux can thus be expressed as,

$$J_V=-D_V\Bigl(\frac{\partial C_V}{\partial x} - \frac{Z^*e E}{kT} C_V\Bigr).$$ (4.19)

The continuity equation for the vacancies gives the following material balance equation,

$$\frac{\partial C_V}{\partial t} + \frac{\partial J_V}{\partial x} + \gamma = 0.$$ (4.20)

where $\gamma$ is a sink/source term which models the recombination or generation of vacancies at sites such as grain boundaries, lattice dislocations, or surfaces. Combining (4.19) and (4.20) we obtain Korhonen's equation,

$$\frac{\partial C_V}{\partial t}-D_V\Bigl(\frac{\partial^2C_V}{\partial x^2}-\frac{Z^*e E}{kT}\frac{\partial C_V}{\partial x}\Bigr)+\gamma=0.$$ (4.21)

This kind of vacancy continuity equation was first considered in [80] for $\gamma=0$ to model the build-up of vacancies at the end of a semi-infinite line.

Instead of an intermediate relationship between the density of lattice sites and hydrostatic stress other researchers [81,82] employed the idea that the vacancy diffusion flux gives rise to volumetric strain which serves to establish stress fields, driving stress migration fluxes. This approach is utilized in [70]. The main advantage is that one does not need the assumption of local stress-vacancy equilibrium, the elastic behavior of metal is not required, and most important, all components of the stress tensor are included. Thus arbitrary geometric shapes in connection with the complex mechanical boundary conditions can be investigated.

4.4.2 Analytical Solution of Korhonen's Equation

Korhonen's equation can be analytically solved for different cases of initial and boundary conditions [51,7,80,83], however from the practical point of view, the most important solution is for the situation where the vacancy flux is blocked at the both ends of a finite line, i.e.

$$\quad J_V(0,t)=J_V(-l,0)=0,$$ (4.22)

on the segment $[-l,0]$, assuming the initial vacancy concentration to be $C_V(x,0)=C_V^{eq}$. In this case, the solution of (4.21) is given in [83] as,

$$\nu(\xi,\zeta)=\frac{C_V(x,t)}{C_V^{eq}}=A_0-\sum_{n=1}^{\infty}A_n (-B_n\zeta+\alpha\xi/2),$$ (4.23)

where,

$$A_0=\frac{C_V(x,\infty)}{C_V^{eq}}=\frac{\alpha \text{exp}(\alpha\xi)}{1-\text{exp}(-\alpha)},$$ (4.24)

is the steady-state solution and,

$$A_n=\frac{16n\pi[1-(-1)^n\text{exp}(\alpha/2)]}{(\alpha^2+4n^2\pi^2)^2}\Bigl(\text{sin}(n\pi\xi)+\frac{2n\pi}{\alpha}\text{cos}(n\pi\xi)\Bigr),$$ (4.25)

$$B_n=n^2\pi^2+\alpha^2/4.$$ (4.26)

For the sake of simplicity, in the equations (4.23)-(4.26) we used the substitutions,

$$\xi=\frac{x}{l},\quad\zeta=\frac{D_Vt}{l^2},\quad\alpha=\frac{Z^*eEl}{kT}.$$ (4.27)

Figure 4.2: Build-up of vacancies at the blocking boundaries for $\zeta =0.01$ (dashed line), $\zeta =0.05$ (points), $\zeta =0.2$ (full line).

Figure 4.3: Build-up of hydrostatic pressure at the blocking boundaries for $\zeta =0.01$ (dashed line), $\zeta =0.05$ (points), $\zeta =0.2$ (full line).

For $\alpha=4$ in Figure 4.2 and Figure 4.3 the distribution of the normalized vacancy concentration $\nu$ and the reduced hydrostatic pressure,

$$\eta=\frac{\Omega\sigma}{kT},$$ (4.28)

are presented for different electromigration stressing times. The stress build-up continues until the electromigration driving force is counter-balanced by intrinsic stress in the line. That is the case in Figure 4.3 for $\zeta =0.2$ where the stress profile is a straight line. The condition for the stress-electromigration equilibrium is,

$$\frac{\Omega}{kT}\frac{\partial\sigma}{\partial x}= \frac{Z^*e\rho J}{kT} .$$ (4.29)

Equation (4.29) is discussed in more detail in [48]. If for the given operating conditions the stress-electromigration equilibrium is reached before a critical stress threshold is build, no damage is produced and the interconnect is virtually ``immortal'' [48,49,47].

4.4.3 The Sink/Source Term $\gamma$

According to the discussion presented in [50], the sink/source term $\gamma$ models the interaction of vacancies between the grain boundary and the grain. The vacancies are annihilated or produced inside the grain boundaries with a rate $\gamma$ which depends on the vacancy concentration in the vicinity of the grain boundaries [50],

$$\gamma = -\frac{C_V-C_V^{eq}}{\tau},$$ (4.30)

where $\tau$ is the average lifetime of the vacancy in the grain boundaries. The effect of the non-zero sink/source term in (4.21) is the prolongation of the time needed for vacancy concentration build-up.

4.4.4 An Explanation of Black's Law

Today a semi-empirical model equation (4.5) is widely used to extrapolate interconnects' time to failure under accelerated test conditions compared to operating conditions. Based on experimental observations [3], the value of the current density exponent $n$ in (4.5) is found to be $n=1$ for the nucleation failure mechanism in which failure is dominated by the time required to build-up a critical stress or a critical vacancy concentration. For the void-growth mechanism, in which failure is determined by the growth of the void to the critical size, the current density exponent is $n=2$.

Considering the nucleation failure mechanism for arbitrary given critical vacancy concentration $C_V^f$, Shatzkes and Lloyd [80] derived following the equation for the time to failure,

$$t_f=2\frac{C_V^f}{D_V^0}\Bigl(\frac{kT}{Z^*e\gamma_T}\Bigr)^2 \frac{1}{J^2} (E_a/kT),$$ (4.31)

with $D_V^0$ the pre-exponential factor for grain-boundary self-diffusivity. The derivation of (4.31) is based on a discussion similar to that given in Section 4.4.1.

4.4.5 Sarychev's Model

Based on the model propositions in [50] Sarychev et al. [70] developed a new model which can be easier applied to realistic interconnect structures. The model connects the evolution of the stress tensor with the diffusion of vacancies whereas the influence of the geometry of the metalization is included. The cause of the stress is the local volume change which is generated by vacancy migration and generation due to electromigration. This mechanical stress opposes electromigration driven vacancy migration. The additional stress load to be included stems from residuum stress as result of the technological process and from mechanical stress due to the thermal mismatch between copper, barrier, and the passivation layer.

The crucial importance for the void-nucleating condition is the dynamics of vacancies in the interconnect metal. The behavior of vacancies can be basically described by the following two equations [7,70],

$$\mathbf{J}_V=-D_V\Bigl(\nabla C_V -\frac{Z^{*}e}{kT\gamma_E}C_V\mathbf{J}+\frac{f\Omega}{kT}C_V\nabla \sigma+\frac{Q^{*}}{kT^2}C_V\nabla T\Bigr),$$ (4.32)

$$\frac{\partial C_V}{\partial t}=-\nabla\cdot\mathbf{J}_V+G(C_V),$$ (4.33)

where $C_V$ is the vacancy concentration, $D_V$ the vacancy diffusivity, $kT$ is the thermal energy, $Z^{*}e$ is the effective valence, $Q^{*}$ is the heat of transport, and $G(C_V)$ is the source function that models vacancy generation and annihilation processes,

$$G(C_V)=L_v \Bigl[kT \frac{C_V}{C_{V,0}}-(1-f)\Omega\sigma\Bigr],$$ (4.34)

where $L_v$ is the rate parameter, $C_{V,0}$ is the vacancy equilibrium concentration in the stress free crystal, and $\sigma=(\sigma_{11}+\sigma_{22}+\sigma_{33})/3$ is the spherical part of the stress tensor.

Local vacancy migration and generation gives rise to the local volume deformation described by the kinetic relation,

$$\frac{\partial \varepsilon^{v}_{kk}}{\partial t}=\Omega[f \nabla \mathbf{J}_V+(1-f)G(C_V)],$$ (4.35)

where $\varepsilon^{v}_{kk}$ denotes the trace of the inelastic strain tensor component due to vacancy generation and migration. Equations (4.32), (4.33) and (4.34) should be solved in conjunction with the electro-thermal equations and the stress equations [70].

$$\mu \nabla^2 u_i + (\lambda + \mu)\frac{\partial}{\partial x_i}(\... ...igl(\frac{\partial u_i}{\partial x_j}+\frac{\partial u_j}{\partial x_i}\Bigr) ,$$ (4.36)

where $\varepsilon_{ij}$ is the total strain tensor, $f$ is the vacancy to atom volume proportion factor, $\lambda$ and $\mu$ are Lame's constants, $B=(3\lambda+2\mu)/3$ is the bulk modulus, and $\mathbf{u}=(u_1,u_2,u_3)$ is the displacement vector.
From the solution of these equations, considering the specific grain structure of the interconnect, we obtain two essential pieces of information. First, we know where in the interconnect geometry void nucleation takes place, secondly we can estimate how much time is needed in the case of a given interconnect geometry and material composition to reach a critical level of tensile stress (or vacancy concentration) at the void nucleation site. The real driving forces of the material transport flux described with (4.32) are the gradients of the electrical field $\nabla\varphi$ and the mechanical stress gradient $\nabla\mathbf{\sigma}$. The mechanical stress gradient is the response of the interconnect metal to the electromigration stress. If the electrical current is so small that it produces an electromigration stress small compared to the mechanical stress, no electromigration takes places. Void nucleation occurs at the moment of equality of these two driving forces.

There has been some difference between the understanding of this sink/source function $G(C_V)$ in the Sarychev model [70,84] and models developed on the basis of Korhonen's work [51,7]. Korhonen [51] assumes that the sink/source function represents the interaction of vacancies between the grain boundary and the grain: annihilation of grain boundary vacancies if their concentration is larger than equilibrium concentration, and their production in the opposite case (Section 4.4.3). On the other hand, Sarychev [70] makes the sink/source term independent of the crystal texture and states that recombination/annihilation takes place everywhere in the bulk.

The described equation system is already solved analytically for special two-dimensional geometry cases in [70]. In [84] a simple finite element solution of the Sarychev's equations for a two-dimensional 20 $\mu$m$\times$ 100 $\mu$m Al thin film is presented.

4.4.6 Finite Element Discretization of the Basic Equation

We are going to demonstrate the numerical handling of the basic continuum equation of the Sarychev model. For the sake of simplicity we assume unpassivated interconnects where the stress relaxation to the hydrostatic state is faster than vacancy diffusion, and where the vacancy concentration in the initial phase of electromigration stressing does not deviate from the equilibrium concentration. In that case the vacancy dynamics is described by a three-dimensional Korhonen-type equation,

$$\frac{\partial C_V}{\partial t}=D_V\Bigl(\Delta C_V+\frac{Z^*e}{kT}\nabla\varphi\nabla C_V\Bigr)-\frac{C_V-C_V^{eq}}{\tau},$$ (4.37)

with the weak formulation,

$$\int_T\frac{\partial C_V}{\partial t}N_j (\mathbf{x}) d\Omega=-D... ...(\mathbf{x}) d\Omega-\int_T \frac{C_V-C_V^{eq}}{\tau}N_j(\mathbf{x}) d\Omega.$$ (4.38)

Figure 4.4: Electrical potential distribution [mV]

Figure 4.5: Vacancy distribution [cm$^{-3}$].

By introducing the spatial discretization of the vacancy concentration $C_V$ and the electrical potential $\varphi$ with the linear basis functions $N_i(\mathbf{x})$,

$$C_V=\sum_{i=1}^4 C_{V,i} N_i(\mathbf{x}),\quad \varphi=\sum_{i=1}^4 \varphi_i N_i(\mathbf{x}),$$ (4.39)

and backward Euler scheme for the time discretization, we obtain,

$$\sum_{i=1}^4\Bigl(\int_T N_i(\mathbf{x}) N_j(\mathbf{x}) d\Omega... ...l(\int_T \nabla N_i(\mathbf{x}) \nabla N_j(\mathbf{x}) d\Omega\Bigr)C_{V,i}^n+$$

$$+D_V\frac{Z^*e}{kT}\sum_{i=1,p=1}^4\Bigl(\int_T \nabla N_i(\mathbf{x}) N_p(\mathbf{x}) \nabla N_j(\mathbf{x}) d\Omega\Bigr)C_{V,i}^n\varphi_{i}^n-$$

$$-\frac{1}{\tau}\Bigl(\sum_{i=1}^4\Bigl(\int_TN_i(\mathbf{x})N_j(\mathbf{x}) d\Omega\Bigr)C_{V,i}^n-C_V^{eq}\int_TN_j(\mathbf{x}) d\Omega\Bigr).$$ (4.40)

The construction of the nucleus matrix follows the scheme described in Section 3.8. The fully discretisatised Sarychev model includes also the finite element discretization of equations (4.6)-(4.11).

Already a simple electromigration model described by (4.37) enables some insight in the electromigration caused material transport phenomena taking place in the copper interconnect layout. In Figure 4.4 the interconnect via with barrier layer is displayed. Here we can also see the electrical potential distribution which determines the electromigration driving force intensity and direction. The vacancy concentration distribution, obtained by solving equation (4.37), shows a peak value at the bottom of the via (Figure 4.5). This peak concentration region indicates a probable location of void nucleation.