Impurity-induced island pinning during electromigration

High current densities in metal conductors can lead to a drift of the metal ions. This effect is known as electromigration (EM). EM is a major reliability issue for interconnects in integrated circuits [1–3] but can also be used to manipulate nanostructures in a desired way [4–6]. Empirically, the reason for the ion drift is an electromigration force $\vec{F}_{\text{EM}}$, which is proportional to the local electric field $\vec{E}$ [7]:

$$\begin{equation} \vec{F}_\text{EM} = Z^*e \vec{E}. \end{equation} \tag{ 1 }$$

In this equation, $Z^*$ denotes the effective valence of the metal ions and e is the elementary charge. The sign of $Z^*$ characterizes the direction of the ion drift. For many metals, like silver (Ag), one finds theoretically and experimentally $Z^*<0$ [7–10]. This implies a drift opposite to the direction of $\vec{E}$.

It has long been known empirically that impurity atoms can reduce electromigration [1,2,11] and therefore increase the lifetime of interconnects. However, EM simulations with more than one atomic species represent a computational challenge, which so far has rarely been addressed, and if so, only with simplified models. For example, in [10,12–16] the electromigration of clusters (islands or voids) is considered only for atoms of one type. Electromigration of an impurity atom in the bulk has been studied with the so-called five-frequency model [17]. By means of ab initio calculations, Bevan et al. [18] showed that certain atomic processes have reduced rates, when impurities are involved. The present paper presents the first material-specific simulation of electromigration of islands on a surface in the presence of an impurity. This is achieved by extending the self-learning kinetic Monte Carlo (SLKMC) method [19–21] for copper (Cu)-contaminated Ag islands on Ag(111).

We generalize our 3D SLKMC model such that the temporal evolution of solids with two kinds of atoms (A and B) can be simulated. In our specific study, atom type A corresponds to Ag and atom type B corresponds to Cu.

In the SLKMC model thermally activated atom hops are performed with a rate ν given by an Arrhenius law:

$$\begin{equation} \nu = \nu_0 \exp{\left(-\frac{E_\text{a}}{k_\text{B}T}\right)}\!. \end{equation} \tag{ 2 }$$

In this equation $k_\text{B}$ is the Boltzmann constant, and T is the absolute temperature. The attempt frequency is assumed to be independent of the atomic process, $\nu_0 = 10^{12}\ \text{s}^{-1}$ [22], and the activation energies $E_\text{a}$ will be calculated from a suitable potential. The attempt frequency and the activation energy have errors. However, as long as the effect on the ratio between any two rates is weak, the error amounts essentially to a rescaling of time. The SLKMC model differs from commonly used KMC implementations [22,23] by combining the accuracy of rates calculated from a realistic potential with the efficiency of a rate catalog. The activation energies are calculated by means of the drag method [24] and saved in a database, from which they are retrieved, when needed during the simulation. In case, a configuration is encountered, which cannot be found in the database yet, the corresponding activation energies are calculated on the fly and added to the rate catalog. To ensure fast addressing within the rate catalog, bit-coding of the hopping processes is used. They are characterized by the occupied lattice sites of two shells surrounding the initial and final position of the hopping atom.

To identify a hopping process within an atom configuration of types A and B, we use five indices: the first one corresponds to the hopping direction, the second and third correspond to the two shells occupied with atoms of type A and the fourth and fifth correspond to the two shells occupied with atoms of type B. An example of this coding is shown in fig. 1. In this example, one atom of type B is embedded within atoms of type A. However, this coding is general enough to consider more than one atom of type B.

Zoom In Zoom Out Reset image size

We implemented an EAM potential [25,26] in order to capture the many-body nature of metallic bonding correctly. We used a well-tested version of it for Ag-Cu-alloys by Williams et al. [27].

In order to simulate EM, we modify the energy term in eq. (2) by $\vec{F}_\text{EM}$ times the vector connecting the initial hopping position and the saddle point $\vec{x}_\text{is}$:

$$\begin{equation} \Delta E_\text{EM} = -\vec{F}_\text{EM} \cdot \vec{x}_\text{is}. \end{equation} \tag{ 3 }$$

The values of $\vec{x}_\text{is}$ are also gathered during the calculation of $E_\text{a}$.

The simulations are performed at $T=350\ \text{K}$. Three parameters are varied. First, the number N of embedded Cu atoms ($N = 0$ or 1), and second, the radius R of the island ($R/r_0 = 3$ to 15, where $r_0 = 0.289\ \text{nm}$ is the nearest-neighbor distance of Ag at room temperature [28]). The third parameter is the ratio of effective valences $Z_\text{Cu}^{*}/Z_\text{Ag}^{*}$. For the physical system considered it has a fixed value, which, however, is not known. Therefore, we consider values in the range $Z_\text{Cu}^{*}/Z_\text{Ag}^{*} \in \left\{-10, -2, 0, 2, 10 \right\}$. Thus, we can also assess qualitatively how the sign and relative strength of the electromigration forces influences the island drift. We also consider negative ratios, because for Ag mass transport in the direction of $\vec{E}$ or $-\vec{E}$ was observed [10,29]. The broad range of ratios reflects the large deviations of experimentally known $Z^*$ [9]. For monolayer islands on a flat metallic surface it is legitimate to assume a constant electric field [14]. Therefore, $Z_\text{Cu}^{*}/Z_\text{Ag}^{*}=| \vec{F}_\text{EM} ^\text{Cu}| / | \vec{F}_\text{EM}^\text{Ag}|$. For Ag atoms we set $|\vec{F}_\text{EM} ^\text{Ag} | = 0.6\ \text{meV}/r_0$ as in [16].

For the investigated systems, the dominant mode of mass transport is edge diffusion. Although considered in the model, an additional terrace diffusion can be neglected [16].

In our simulations we are interested in the drift velocity parallel to $\vec{F}_\text{EM}$. As we set $\vec{F}_\text{EM} = |\vec{F}_\text{EM} | \cdot \vec{e}_x$ (where $\vec{e}_x$ is the x-unit vector) we calculate the x-component of the drift velocity v_x. To do so, the temporal evolution of the x-center-of-mass coordinate r_x of the islands is investigated.

Every simulation of each set of parameters is repeated ten times. Subsequently, the ten $r_x(t)$ trajectories are averaged to calculate the mean velocity $\langle v_x\rangle$.

The mean drift velocities and the corresponding error bars are shown in fig. 2 for each set of parameters. The results of pure Ag islands are in accordance to the ones found in [16]. For larger islands $(R/r_0 > 15)$ we expect $\langle v_x\rangle\propto 1/R$ [16]. If an Ag island contains one Cu atom, our simulations reveal that $\langle v_x\rangle$ decreases. The ratio $v_x^\text{Ag}/v_x^\text{Ag-Cu}$ depends on R and becomes larger with decreasing R: At $R/r_0 = 3$ and $Z_\text{Cu}^{*}/Z_\text{Ag}^{*}=-10$ a maximal ratio of approximately 40 is observed. The drift reduction is only weakly influenced by $Z_\text{Cu}^{*}/ Z_\text{Ag}^{*}$: qualitatively, higher values lead to higher drift velocities. The results for $R/r_0 = 4$ seem to be outliers. But as Heinonen et al. found in [30], the diffusion constant D shows size-dependent oscillations for small islands. Since $v_x \propto D$, we see such an oscillation.

Zoom In Zoom Out Reset image size

The reason for the drift reduction induced by a Cu impurity can be understood by investigating single $r_x(t)$ trajectories. Two representative trajectories are shown in fig. 3 ($R/r_0 = 5$, $Z_\text{Cu}^{*}/Z_\text{Ag}^{*} = 10$ or −10). In addition, the trajectory of a pure Ag island without a Cu impurity is plotted.

Zoom In Zoom Out Reset image size

In contrast to the pure Ag island, the trajectories of both Ag-Cu islands show plateaus. This is due to a pinning of the whole island.

The pinning lasts for a considerable but finite time (plateaus). After this pinning time τ, a depinning occurs. Once the depinning happened, both Ag-Cu islands migrate with a mean drift velocity $\langle v_x^\text{depin}\rangle= 14.75~r_0/\text{s} \pm 1.13~r_0/\text{s}$ which is within the standard error of the mean drift velocity of pure Ag islands $\langle v_x \rangle= 14.22~r_0/\text{s} \pm 0.50~r_0/\text{s}$. The pinning time τ and migration range between two plateaus $\Delta r$ depend on $Z_\text{Cu}^{*}/ Z_\text{Ag}^{*}$. The plotted trajectories show an increase of τ and a decrease of $\Delta r$ with $Z_\text{Cu}^{*}/Z_\text{Ag}^{*}$, respectively. Therefore, we see a dependence of $\langle v_x\rangle$ on the ratios as is shown in fig. 2.

To understand the pinning effect in more detail, we discuss five representative simulation snapshots of an Ag-Cu island with $R/r_0 = 5$ and their qualitative temporal occurrence (see fig. 4).

Zoom In Zoom Out Reset image size

A simulation starts at a time t₀ with a round island as shown in fig. 4(a). The Cu atom is placed at an arbitrary site within the island.

$\vec{F}_\text{EM}^\text{Ag}$ leads to a drift of Ag atoms to the right via edge diffusion. Therefore, the whole island migrates rightwards until the Cu atom becomes the second nearest neighbor of a left-edge atom. This case is shown in fig. 4(b) and occurs at a time t₁. The difference between t₀ and t₁ is rather short and is determined by $\langle v_x^\text{Ag}\rangle$. As mentioned before, our SLKMC model considers two shells surrounding a hopping atom. Therefore, the Cu atom does not influence the edge diffusion until t₁. After t₁, the Cu atom influences the hopping processes of its surrounding Ag atoms. It turns out that many calculated activation energies for configurations including a Cu atom are higher than for a pure Ag configuration (note that for different impurities this may not be the case [18]). Therefore, the pinning starts at t₁.

Due to random fluctuations the Cu atom can become an edge atom. This case is shown in fig. 4(c) and occurs at a time t₂. Then the Cu atom can diffuse along the island edge which leads to depinning. Whether and how far it diffuses, depends on the ratio $Z_\text{Cu}^{*}/Z_\text{Ag}^{*}$. In our simulations, the Cu atom is strongly bound to the island, so that even for $\left|Z_\text{Cu}^{*}/Z_\text{Ag}^{*}\right|=10$ we never observed that it detaches from the island edge. Thus, a variation of $Z_\text{Cu}^{*}/ Z_\text{Ag}^{*}$ only affects the distance, the Cu atom diffuses, before it becomes surrounded again by Ag atoms. The larger $Z_\text{Cu}^{*}/Z_\text{Ag}^{*}$, the more likely a Cu edge atom will hop in the direction of $\vec{F}_\text{EM}^\text{Ag}$. For $Z_\text{Cu}^{*}/Z_\text{Ag}^{*} \gg 1$ it is very likely that the Cu atom reaches the opposite island edge. Thus, we see $\Delta r \approx 2 R = 10 r_0$ for $Z_\text{Cu}^{*}/Z_\text{Ag}^{*} = 10$ in fig. 3. In the case of $Z_\text{Cu}^{*}/Z_\text{Ag}^{*} \ll 1$ the migration range $\Delta r$ becomes small or even zero. We see $\Delta r = 0$ if the Cu atom becomes surrounded again, before it diffuses in the direction of $\vec{F}_\text{EM}^\text{Ag}$.

In fig. 4(d) the Cu atom had diffused about $4/3 R$ along the island edge. The Cu diffusion ends at a time t₃. The difference $t_3-t_2$ is small compared to $t_2-t_1$ because edge diffusion occurs at short time scales. The pinning time is characterized by $\tau = t_3 - t_1$.

After t₃, when the Cu atom becomes embedded into the island again, the Ag atoms migrate freely with $\langle v_x^\text{depin}\rangle\approx\langle v_x^\text{Ag}\rangle$. The snapshot in fig. 4(e) shows the Ag diffusion after the depinning at time t₄. The pure Ag diffusion will last until the Cu atom becomes the second nearest neighbor of a left-edge atom. If this happens, a new pinning occurs and the described procedure takes place again.

We have discussed the case $R/r_0 = 5$. For other choices of R, we can also see a pinning and depinning of the island. But with increasing R it becomes more difficult to separate between the pinning and pure Ag drift since $\langle v_x^\text{Ag} \rangle$ decreases [16]. The r_x trajectories approximate the ones of pure Ag islands and the difference between $\langle v_x^\text{Ag} \rangle$ and $\langle v_x^\text{Ag-Cu} \rangle$ becomes smaller as is shown in fig. 2.

Using an extended SLKMC model, we were able to simulate the electromigration-induced drift of monolayer Ag islands on Ag(111) which contain one Cu atom. Our simulations show that the single Cu atom causes a pinning of the whole island. The pinning occurs regardless of the value and sign of $Z^{*}_\text{Cu}/Z^{*}_\text{Ag}$ and induces a reduction of the island drift velocity. The difference between $v_x^{\text{Ag}}$ and $v_x^{\text{Ag-Cu}}$ is very large for small R but decreases, if R grows. This effect is also expected for different alloys and island structures. The pinning described in this paper depends crucially on the increase of activation energies due to an impurity nearby, such as already calculated by Bevan et al. [18] in a different context. For which alloys the effect will occur can be predicted by density functional theory (DFT) calculations of these activation energies.

For future work it would be interesting to investigate Ag islands with more than one Cu impurity. This is possible with our extended SLKMC model, since it is general enough to consider more than just one Cu atom and is only limited by the computation time to calculate new activation energies. In addition, it is possible to change the used potential and investigate the drift of different alloys. Another direction, in which progress seems possible, is to generalize existing phenomenological theories [14,15] in order to include the position of the impurity relative to the island edge and its influence on the temporal evolution.

Helpful discussions with L. Brendel are gratefully acknowledged.