Cu(001) surface alloy

Cu(001) surface alloy

Diffusion-mediated processes in Pt/Cu(001) surface alloy Journal Pre-proof Diffusion-mediated processes in Pt/Cu(001) surface alloy S.A. Dokukin, S...

9MB Sizes 0 Downloads 23 Views

Diffusion-mediated processes in Pt/Cu(001) surface alloy

Journal Pre-proof

Diffusion-mediated processes in Pt/Cu(001) surface alloy S.A. Dokukin, S.V. Kolesnikov, A.M. Saletsky, A.L. Klavsyuk PII: DOI: Reference:

S0039-6028(19)30671-5 https://doi.org/10.1016/j.susc.2019.121515 SUSC 121515

To appear in:

Surface Science

Received date: Revised date: Accepted date:

4 September 2019 23 October 2019 24 October 2019

Please cite this article as: S.A. Dokukin, S.V. Kolesnikov, A.M. Saletsky, A.L. Klavsyuk, Diffusion-mediated processes in Pt/Cu(001) surface alloy, Surface Science (2019), doi: https://doi.org/10.1016/j.susc.2019.121515

This is a PDF file of an article that has undergone enhancements after acceptance, such as the addition of a cover page and metadata, and formatting for readability, but it is not yet the definitive version of record. This version will undergo additional copyediting, typesetting and review before it is published in its final form, but we are providing this version to give early visibility of the article. Please note that, during the production process, errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain. © 2019 Elsevier B.V. All rights reserved.

Highlights • Order-disorder phase transition in Pt/Cu(001) surface alloy is observed at 350-400 K • Dissolution of small Pt clusters in the topmost layer of Cu(001) is investigared • Suppression of vacancy clusters electromigration in Pt/Cu(001) alloy is predicted • Universal TB-SMA potential parameters for Pt/Cu(001) and Pt/Cu(111) are fitted

1

2

Diffusion-mediated processes in Pt/Cu(001) surface alloy S.A. Dokukina,b,∗, S.V. Kolesnikova , A.M. Saletskya , A.L. Klavsyuka a

Faculty of Physics, Lomonosov Moscow State University, Moscow 119991, Russian Federation b A.M. Obukhov Institute of Atmospheric Physics, Russian Academy of Sciences, Pyzhevsky per., 3, 119017, Moscow, Russia

Abstract In this paper we present the investigation of the diffusion-mediated processes in a Pt/Cu(001) surface alloy. We use the interatomic semiempirical TBSMA interatomic potentials and the self-learning kinetic Monte Carlo method to investigate the following processes: order-disorder phase transition in the Pt/Cu(001) surface alloy, dissolution of small Pt clusters, and electromigration of small vacancy clusters in the topmost layer of the Pt/Cu(001) surface alloy. We have fitted the universal parameters of the interatomic TB-SMA potentials for the Pt–Cu system. The potentials reproduce bulk properties of copper and platinum as well as properties of the Pt/Cu(001) and Pt/Cu(111) surface alloys. Keywords: metals and alloys, nanostructured materials, transition metal alloys and compounds, atomic scale structure, kinetics, computer simulations 1. Introduction Surface alloys are intensively used for catalysis [1], improvement of the mechanical properties of surfaces [2], construction of fuel cells [3] etc. Pt/Cu(001) surface alloy formation was investigated in experimental [4–10] and theoretical [11–14] studies. Several important features were noticed in ∗

Corresponding author Email address: [email protected] (S.A. Dokukin)

Preprint submitted to Elsevier

October 24, 2019

investigations of the Pt/Cu(001) surface alloy formation with high concentration of Pt atoms. Formation of the c(2 × 2) structure is observed after deposition of 0.5 – 2 ML of Pt atoms on a Cu(001) substrate [4–6]. Annealing of this alloy results in formation of the c(2×2) layer capped with the pure Cu p(1 × 1) layer [8–11]. Growth of a pure Pt thin film is observed at higher concentrations of Pt atoms [5]. However, experimental studies of the Pt/Cu(001) surface alloy formation with low concentration of Pt atoms are not presented in the literature. Several atomistic properties of the Pt/Cu(001) surface alloy were revealed in theoretical studies. Calculations with the semiempirical Bozzolo-Ferrante-Smith potentials showed that the energy of the Pt dimer in position of the third neighbors (p(2×2) structure) is lower than in position of the second neighbors (c(2 × 2) structure) [12]. Ab initio calculations revealed the most energetically favourable positions, density of states and absorption spectra of Pt atoms on the Cu(001) surface [13, 14]. Two deposition modes are used in experimental studies. First is the deposition of single atoms. Second is the deposition of small clusters [15, 16]. Such clusters can embed into substrate upon further annealing [17]. The clusters begin to gradually dissolve if the mixing enthalpy of alloy is negative. Dissolution of clusters was observed in such systems as B clusters in Si [18] and Cu nanoparticles in SiO2 [19]. The kinetic Monte Carlo (kMC) method can be successfully applied to the simulation of the dissolution process [20]. Another important problem is investigation of the electromigration processes. Electromigration is the source of destructive processes in metals [21, 22]. However, incorporation of impurity atoms can significantly reduce the rate of electromigration [23]. The kMC method can be used to calculate drift velocity of islands caused by electromigration [24–28]. Using the self-learning kMC (SLkMC) algorithm allows to increase the precision of calculations and consider multi-component systems. The SLkMC simulations revealed the nonmonotonic dependence of the drift velocity on temperature and size of islands in Ag(111) and Ag(001) films [25], reduction of the island drift velocity after introduction of a single Cu atom into an Ag island on the Ag(111) surface [26] and oscillations in the dependence of the drift velocity on the size of a vacancy clusters in the Cu(001) surface [27, 28]. A lot of theoretical methods are used to investigate the physical parameters of surface alloys. Ab initio methods [29–32] are the most reliable but also quite complicated and time-consuming. Classical approach with semiempirical potentials allows to significantly accelerate calculations. Potentials based on the tight binding method in the second moment approximation (TB4

SMA) [33, 34] can simultaneously reproduce the atomic adsorption energies, the values of the diffusion barriers and precisely describe the bulk structure of substrates. The parameters of the potentials can be fitted using the values of the binding energies and the diffusion barriers to obtain agreement between the experimental results and the simulations [35–37]. The kMC method [38] allows to simulate the formation of surface alloys at experimental time scale. Application of the SLkMC [39–43] technique, in which diffusion barriers are calculated on-the-fly, significantly increases accuracy of the calculations and allows to consider multi-component alloys. Reliability of this approach was shown for various systems: Cu clusters on Cu(111) [39], Al atoms in Mg [40], Ag clusters on Ag(111) [41], vacancies in Si [42], Co nanostructures in Cu(001) [43]. The formation and the physical properties of the Pt/Cu(001) surface alloy with low concentration of Pt atoms is poorly studied. We employ the SLkMC technique and the TB-SMA interatomic potentials to the investigation of properties of the Pt-Cu surface alloy. The following interesting problems are discussed in the paper: (i) the formation of the Pt/Cu(001) surface alloy under various conditions, (ii) the dissolution of the Pt cluster in the topmost layer of the Cu(001) substrate, and (iii) the influence of Pt atoms on the electromigration rate of the vacancy cluster in the topmost layer of the Cu(001) substrate. The paper is organized as follows. The computational method is discussed in Section 2. We present and discuss the simulation results in Section 3. Section 4 concludes our paper. The details of the fitting procedure and some other technical details are presented in Supplementary Materials. 2. Computational method We use the TB-SMA potentials [33–35] to simulate the interatomic interactions. The attractive term Ebi (band energy) contains the many-body interactions. Pair interactions (Born-Mayer form) describe the repulsive part Eri . The cohesive energy EC is the sum of the band energy and the repulsive part: X EC = (Ebi + Eri ), (1) i

v    uX u rij i 2 ξαβ exp −2qαβ − 1 fc (rij ), Eb = −t r0αβ j 5

(2)

Eri

      X rij rij 0 1 − 1 + Aαβ exp −pαβ − 1 fc (rij ), = Aαβ r0αβ r0αβ j

(3)

where rij is the distance between the atoms i and j; α and β are types of the atoms; r0αβ , A0αβ and A1αβ are adjustable parameters of an interatomic interaction; pαβ and qαβ describe the decay of the interaction strength with the distance between the atoms; ξαβ is an effective hopping integral. The cut-off function fc (rij ) has the form [44]  1 , rij ≤ Ron   2 2 2 )2 (R2 2 (Rof f −rij of f −3Ron +2rij ) , Ron < rij < Rof f fc (rij ) = (4) 2 2 3 (Rof f −Ron )   0 , rij ≥ Rof f where Ron = 6.0 ˚ A and Rof f = 6.5 ˚ A are cut-off distances. Parameters for the Cu-Cu interaction were taken from Ref. [45]. These potentials were previously used in multiple investigations of the formation of various atomic structures (supported islands, surface and bulk vacancies, nanocontacts, dendrites) [46–50]. The experimental data and the ab initio calculations were used to fit the parameters for the Pt-Pt and Pt-Cu interactions. The fitting procedure is similar to the one described in Ref. [37]. Details of the fitting procedure are discussed in the Supplementary Materials. The parameters of the TB-SMA potentials are presented in Table 1. Table 1: Parameters of interatomic interactions.

Parameter A1 (eV) A0 (eV) ξ (eV) p q r0 (˚ A)

Cu-Cu[45] 0.000 0.085 1.224 10.939 2.280 2.556

Cu-Pt 0.240 0.230 2.033 10.597 3.606 2.716

Pt-Pt -0.216 0.107 1.904 10.425 4.116 2.986

The SLkMC method [51] was applied for the realistic simulation of the Pt/Cu(001) surface alloy formation. The SLkMC calculation starts with an empty database and only necessary diffusion barriers are computed. In our method the nudged elastic band (NEB) [52] method is integrated in the 6

SLkMC for the calculation of the diffusion barriers. The positions of the Pt and Cu atoms are determined in a fully relaxed geometry. The slab used to calculate the barriers consists of eight layers with 200 atoms per layer 2 (36 × 36 ˚ A ). Positions of atoms from the two bottom layers are fixed, and periodic boundary conditions are applied to the surface plane. The following formula is used in the kMC algorithm to calculate the rates of the events   Ed 0 , (5) νi = νi exp − kB T where νi0 is the frequency prefactor, Ed is the diffusion barrier, kB is the Boltzman constant and T is the temperature of the substrate. The frequency prefactors for diffusion of single atoms and dimers are equal to 1.5 · 1013 s−1 and 3·1014 s−1 , respectively [51]. The computational cell is a one-layer square 2 grid with 100 × 100 atoms (256 × 256 ˚ A ). Periodic boundary conditions are applied to the surface plane. Maximum concentration of the Pt atoms is 0.5 ML. The concentration of vacancies is much lower than the concentration of atoms in the simulation. Thus, we consider vacancy jumps instead of atomic jumps to accelerate the calculations. We use the random number generator from [53] to improve the accuracy of the calculations. We summarize some data from fitting and validation of parameters of the TB-SMA potentials in Table 2 for the convenience of readers. 3. Results and Discussions 3.1. Order-disorder phase transition In this section we consider the order-disorder phase transition in the Pt/Cu(001) surface alloy at 300–400 K. First, let us briefly discuss the formation of a PtCu alloy in the first layer of the Cu(001) surface. When a Pt atom is deposited to the Cu(001) surface, it can embed into the surface via exchange with Cu atom from the substrate. The diffusion barrier of this event is 0.527 eV (Ed3 in Table 2). The rate of this event is 4.2 · 105 s−1 (6.7 · 107 s−1 ) at 300 K (400 K). As we will show later, this rate is much higher than the characteristic relaxation rates of the alloy. The barrier for the reverse event is relatively high: 1.048 eV (Ed4 in Table 2). Therefore, the rate of such events is low ν = 7.5 · 10−4 s−1 (19 s−1 ) at 300 K (400 K). These events are very rare and can be ignored at considered temperatures. The embedded Pt atom can move in the topmost layer of the substrate via 7

Table 2: Characteristic energies and diffusion barriers discussed in the text of the article. Ed2 is the diffusion barrier for the exchange of the Pt atom with vacancy inside the topmost layer of the Cu(001) substrate. Ed3 and Ed4 are the diffusion barriers for the direct and inverse exchanges of the Pt adatom on the Cu(001) substrate and Cu atom from the (001) (001) substrate. E2in1 and E2in2 are the binding energies of the heterogeneous Pt dimers in the topmost layer of the Cu(001) substrate in the first and the second nearest neighbour position, respectively. See more details in Supplementary Materials.

Quantity Ed2 Ed3 Ed4 (001) E2in1 (001) E2in2

DFT data 0.535 eV 0.527 eV 1.048 eV 0.156 eV 0.038 eV

TB-SMA value 0.620 eV 0.614 eV 1.368 eV 0.326 eV 0.179 eV

exchanges with vacancies. The value of the diffusion barrier of this event is 0.535 eV (Ed2 in Table 2). The rates of such events depend on the concentration of vacancies in the topmost layer of the Cu(001) substrate. The concentration of vacancies in the topmost layer of clean Cu(001) surface can be estimated by the formula [48]   GF , (6) nvac = exp − kB T where kB is the Boltzmann constant, T is the substrate temperature, and GF = 0.393 eV is the Gibbs free energy of the vacancy formation in the topmost layer of the substrate. A detailed analysis of the initial stage of the Pt/Cu(001) alloy is presented in Ref. [12]. However, the characteristic ordering time in the alloy is several orders of magnitude longer than the time of the embedding of a Pt atom. Therefore, this analysis does not allow to fully understand the processes of alloy formation. We use the following simple SLkMC model to study the order-disorder phase transition in the Pt/Cu(001) alloy. We assume that (i) Pt atoms are randomly distributed in the topmost layer of the Cu(001) substrate at the initial moment of time (t = 0), (ii) Cu atoms displaced on the surface do not affect the diffusion of platinum atoms, (iii) diffusion of the 8

embedded Pt atoms occurs due to the exchanges with vacancies. There is only one vacancy in the calculation cell in order to accelerate the simulations. Therefore, the concentration of vacancies in the simulations is n ˜ vac = 10−4 ML (size of the calculation cell is 100×100 atoms). The value of n ˜ vac is higher than the concentration of vacancies nvac calculated with formula (6) at 300– 400 K. Thus, the simulated time t˜ should be rescaled as t = t˜ · n ˜ vac /nvac . The DFT calculations show that Pt atoms in the second nearest neighbour position repel each other (see Table 2). Two Pt atoms in the third nearest neighbour position repel weaker than two Pt atoms in the second nearest neighbour position1 . At first glance, this ratio should result in the formation of the p(2 × 2) surface alloy when the concentration of Pt atoms is equal to 0.25 ML. However, the repulsion energy of a Pt dimer in the second nearest neighbour position is not very high. Therefore, it may be energetically favourable for high concentration of Pt atoms to form locally structured c(2×2) surface alloy due to the many-body nature of the TB-SMA potentials and the relaxation effects. We calculated the Pt-Pt correlation functions at 300 K, 350 K and 400 K and the concentration of Pt atoms 0.25 ML to define the type of the atomic ordering in the Pt/Cu(001) surface alloy. We define the Pt-Pt correlation function as [54] * + X 1 ρ(r0 )ρ(r0 + r) − nP t , (7) C(r) = NM L 0 r where the density ρ(r) is defined to be 1 if the site is occupied by Pt atom and 0 in the other cases, nP t is the concentration of Pt atoms, and NM L is the number of atoms in the monolayer, the angular brackets represent the averaging over all directions of r. We assume that the Pt atoms are randomly distributed after the deposition and, therefore, C(r, t = 0) = 0. The type of the atomic ordering in the alloy at t > 0 is characterized by the extrema points of the correlation function. Fig. 1 shows the correlation functions at 300 K averaged over 100 simulations. The sign of the correlation function for 10 nearest neighbours can be precisely defined at 40 s after the deposition. The values of the correlation function for 50 nearest neighbours (Fig. 1 shows 20 nearest neighbours) are clearly observed at 4000 s after the 1

The binding energies of the Pt dimers in the Cu(001) surface are equal to E2in1 = 0.326 eV, E2in2 = 0.179 eV (see Table 2) and E2in3 = 0.076 eV according to the calculations with the TB-SMA potentials.

9

Figure 1: The Pt-Pt correlation function at 40 s, 4000 s and 40000 s after the deposition of platinum atoms. Concentration of Pt atoms is 0.25 ML. Temperature of the substrate is 300 K. The dashed lines show the Pt-Pt correlation functions in the case of ideal c(2 × 2) (divided by factor 3) and p(2 × 2) (divided by factor 2.5) structures.

10

deposition. The situation at 40000 s is similar to the situation at 4000 s. Therefore, we can suppose that the distribution of atoms in the Pt/Cu(001) alloy reaches the equilibrium at t = 4000 s. The correlation functions for the c(2 × 2) and the p(2 × 2) atomic ordering in ideal (001) surface are shown in Fig. 1 with the dashed lines for comparison. Values of these correlation functions are divided by factors 3 (2.5) for the c(2 × 2) (p(2 × 2)) structures. It is clearly seen from Fig. 1 that the sign of the correlation function of the first 20 nearest neighbours is the same as for the correlation function of the ideal c(2 × 2) structure. It means that Pt atoms arrange in clusters with the c(2×2) structure. At the same time, the number of atoms in the third nearest neighbour position is higher than in the second nearest neighbour position. We conclude that the alloy is composed of the areas with the c(2×2) structure separated by the the small areas with the p(2 × 2) structure. We can approximate the absolute value of the correlation function with the following function [54] in order to check the presence of the long-range order in the alloy: A (8) |C(r)| = δ e−r/ξ + C∞ , r where ξ is the correlation radius, A = const and C∞ = lim |C(r)|. The conr→∞ stant value C∞ plays the role of the parameter of order of the order-disorder phase transition. We can not observe values of the correlation function at r → ∞ due to the limitations of the calculation cell. Moreover, parameter δ is much lesser than 1 for 2D systems [55]. It is possible to approximate the function 1 (9) ln (|C(r)| − C∞ ) = ln A − ln r, ξ with the linear functions in range [rmin , rmax ] which is chosen in a way to exclude fluctuations |C(r)| at short and long distances. The value of C∞ can be either equal to 0 or higher than 0. If the approximation (9) with C∞ = 0 is better, then there is no long-range order in the alloy. Otherwise, C∞ 6= 0 and there is a long-range order in the alloy. There is no long-range order at all temperatures in the range 300-400 K, when the concentration of Pt atoms is 0.25 ML. However, the long-range order is observed for the concentrations of Pt atoms higher than 0.3 ML (see Fig. 2). The correlation radius ξ increases rapidly before and decreases rapidly after the long-range order is established in the alloy in the phase transition region. 11

Figure 2: The time dependencies of the correlation radius 350 K and 400 K, and the concentration of Pt atoms of 0.35 ML and 0.5 ML.

12

Figure 3: The relaxation time as the function of the Pt concentration at 350 K and 400 K. The solid lines show the linear approximations.

13

At the beginning of the simulation (t = 0) the alloy is in a disordered state (C(r) = 0). We define the relaxation time τ as the time of establishment of the long-range order in the alloy. An increase in the concentration of Pt atoms results in the increase of the relaxation time τ , while an increase in the substrate temperature results in the decrease of the relaxation time τ (Fig. 3). The long-range order was not observed during simulations at room temperature. This fact gives the lower limit of the relaxation time at 300 K: τ & 105 s. Fig. 3 shows the relaxation time at 350 and 400 K. The relaxation time error ∆τ can be defined as the peak width at half-height of the ξ(t). The function ln τ (nP t ) can be approximated with a linear function ln τ = a · nP t + b, where a = 0.136 ML−1 and b = 1.682 at T = 350 K, and a = 0.207 ML−1 and b = −5.069 at T = 400 K. The relaxation times τ can be much longer if we take into account the interaction of platinum atoms with copper atoms displaced on the surface. Therefore, our model gives lower estimates of the relaxation times. 3.2. Dissolution of small Pt clusters Let us discuss another mode in which atomic clusters are deposited to the substrate [15, 16]. These clusters can embed into the Cu(001) substrate [17]. The Pt clusters dissolve in the Cu substrate, because it is energetically unfavourable for Pt atoms to be located in the first nearest neighbour position in a Cu substrate. Let us discuss the dissolution of small Pt clusters in the Cu(001) surface using the following simple model. We assume that (i) a round single-layer Pt cluster with a radius of R0 is located in the topmost layer of the Cu substrate in the initial moment of time (t = 0); (ii) this cluster has the fcc structure with a lattice constant equal to the lattice constant of Cu (3.615 ˚ A), and platinum atoms are located in the sites of the Cu substrate; (iii) Cu atoms displaced in the upper layer do not affect the diffusion of platinum atoms; (iv) the diffusion of the embedded atoms occurs via exchanges with surface vacancies. We consider the platinum clusters with initial radii R0 from 15 ˚ A to 40 ˚ A in the simulations. Such radii correspond to the average concentration of Pt atoms from 1% to 8% in the calculation cell2 . Formation of the disordered alloy is observed at such concentrations of platinum atoms 2

We use the one-layer SLkMC model described in the Section 2. Thus, the volume concentration of atoms is equal to the surface concentration of atoms in the framework of our model.

14

in the temperature range of 300–400 K and t → ∞ (see Sec. 3.1). The distribution of atoms in such alloy is close to the random distribution. We perform the simulation with one vacancy in the calculation cell and rescale the simulation time according to the formula (6) similarly to the previous Section.

Figure 4: The radial distribution of the Pt atoms in the Cu(001) surface during the dissolution of the Pt cluster. The initial radius of the cluster is 40 ˚ A. The temperature of the substrate is 350 K. The colored points are the calculated distributions averaged over 100 simulations. The solid lines are the ideal distributions of Pt atoms in the p(1 × 1), c(2 × 2), p(2 × 2), and the fully disordered Pt/Cu(001) alloys (see also Supplementary Materials).

Figure 4 shows the radial distribution Nk (rk ) of the Pt atoms in the first layer of the Cu(001) surface during the dissolution of the Pt cluster, where Nk is the number of Pt atoms in a ring of radius rk = kr0 and width r0 √ centred in the center of the cluster, r0 = a/ 2 = 2.556 ˚ A is the nearest neighbour distance and k ∈ N. The initial radius of the cluster is R0 = 40 ˚ A. 15

The average concentration of the Pt atoms in the calculation cell is 7.72%. The points in Fig. 4 show the distributions obtained by the simulations of the cluster dissolution at 350 K. Each point is obtained by averaging over 100 simulations. Errors of calculations are a few per cent of a value and are not shown in Fig. 4. The solid lines in Fig. 4 show the radial distributions in the ideal p(1 × 1), c(2 × 2), p(2 × 2) alloys and in the random distribution of the Pt atoms with a concentration of 7.72%. The cluster has a “dense core” in which Pt atoms are located at the distances of the first or the second nearest neighbours at the initial stage of the dissolution (t < 50 s). The points between the straight lines p(1 × 1) and c(2 × 2) in Fig. 4 correspond to this region. We define the boundary of the region with a dense core R2 as the intersection of the radial distribution of the Pt atoms and the c(2 × 2) line. Pt atoms are located at the distances of the second or the third nearest neighbours at r > R2 from the center of the cluster. We will call this region a “loose structure”. The points that lie between the straight lines c(2 × 2) and p(2 × 2) in Fig. 4 correspond to this region. We define the outer boundary of the region with a loose structure R3 as the intersection of the radial distribution of the Pt atoms and the p(2 × 2) line. The cluster has a disordered structure at distances r > R3 from its center. We define the outer boundary of the region with a disordered structure Rrand as the intersection of the radial distribution of the Pt atoms and the line corresponding to the random distribution of the Pt atoms. We assume that the cluster boundary is determined by the value of Rrand . The dissolution of cluster occurs as follows. The dense core initially expands, after that it decreases in size, and, finally, disappears at the moment of time t = τ2 (R2 (τ2 ) = 0). The cluster grows in size during t ∈ [0, τ2 ]. This growth is, obviously, the result of the strong repulsion between the Pt atoms (see Table 2). The center of the cluster has a loose structure described in the previous section (i.e. distribution of Pt atoms is the same as in the c(2 × 2) alloy) at t > τ2 . The loose structure of the cluster disappears at the moment t = τ3 (R3 (τ3 ) = 0). The time τ3 is an order of magnitude larger than the time τ2 because the repulsion between the Pt atoms in the second neighbours is lower than in the first neighbours. The cluster has a disordered structure and behaves like a two-dimensional cloud of Brownian particles at t > τ3 . The cluster completely dissolves (Rrand → ∞) after t ∼ 1 h. The functions R2,3 (t) are nonmonotonic and tend to zero after the initial growth of the cluster. At the same time, the function Rrand (t) increases monotonically. 16

Figure 5: The time dependencies of the radii R2 and R3 of the Pt cluster at 350 K. The initial radius of the cluster is 40 ˚ A. The colored points are the calculated distributions averaged over 100 simulations. The solid lines are the approximations with the linear functions R2,3 = a2,3 t + b2,3 . The upper inset shows the dependencies of times τ2 and τ3 on the inverse temperature of the cooper substrate. The solid lines in the upper inset 0 are the approximations with functions τ2,3 = τ2,3 exp(E2,3 /kB T ). Lower inset shows the dependencies of times τ2 and τ3 on the initial radius of the cluster at 300 K. Solid lines in 0 the downer inset are the approximations with the functions τ2,3 = τ¯2,3 exp(α2,3 R0 ).

17

Let us discuss in more detail the dynamics of the cluster core after the initial growth. Fig. 5 shows the dependencies of the radii of the dense core R2 and the loose region R3 on the time at different temperatures of the copper substrate in the range of 300–350 K. The initial cluster radius is R0 = 40 ˚ A. The points show the averaged values of the radii. The solid lines are the linear functions R2,3 = a2,3 t + b2,3 , where the coefficients a2,3 and b2,3 are obtained with the least squares method. It is seen from Fig. 5 that the functions R2,3 (t) decrease proportionally to the time3 . The characteristic times of the cluster dissolution can be easily calculated with the formula τ2,3 = −b2,3 /a2,3 if the values of the coefficients a2,3 and b2,3 are known. The upper inset in Fig. 5 shows the dependencies of the τ2 and τ3 on the inverse temperature. It can be seen that these dependencies in the temperature range 300-350 0 exp(E2,3 /kB T ), where K can be approximated with the functions τ2,3 = τ2,3 −14 0 s, and E3 = 1.115 eV, τ30 = 3.6 × 10−14 s. E2 = 1.134 eV, τ2 = 1.0 × 10 The characteristic times τ2 and τ3 depend exponentially on the initial radius of the cluster R0 . The lower inset in Fig. 5 shows the dependencies of τ2 and τ3 on the radius R0 in the range from 15 ˚ A to 40 ˚ A at 300 K. It can be seen that these dependences can be approximated with the functions τ2,3 = −1 −1 0 A , exp(α2,3 R0 ), where α2 = 0.135 ˚ A , τ¯20 = 5.3×102 s, and α2 = 0.128 ˚ τ¯2,3 τ¯20 = 1.9 × 103 s. Figure 6 shows that the Pt atoms located on the periphery of the Pt cluster behave like a cloud of non-interacting Brownian particles. The mean square displacement of a Brownian particle is proportional to time hr2 i = 4Dt in the two-dimensional case, where D is the diffusion coefficient. Therefore, the radius√of the Brownian particle cloud √ increases as a square root of time Rrand ( t) can be approximated with the Rrand ∼ t. The dependencies √ functions Rrand = a ˜ t + ˜b in the temperature range of 300–350 K. The insert in Fig. 6 shows the dependence of the diffusion coefficient D ≈ a ˜2 /4 on the inverse temperature 1/kB T for the cluster with the initial radius of R0 = 40 ˚ A. It is clearly seen that this dependence can be approximated 3

The origin of this linear dependence can be easily understood in the framework of the simplest analytical dissolution model. If N is the number of atoms in the 2D cluster, then √ the rate of the dissolution is proportional to the perimeter of the cluster P ∼ N . The time dependence√of the number of atoms in the cluster can be found from the equation dN (t)/dt = −α N , α = const with the initial condition N (0) = N0 . The solution of 2 this equation is N (t) = N0 − αt2 /2 . Therefore, the radius of the cluster is R(t) = p √ r0 N (t)/π and, finally, R(t) = R0 − t · αr0 /2 π.

18

Figure 6: The radius Rrand of the Pt cluster as the function of the square root of time at 300–350 K. The colored points are the calculated distributions averaged√over 100 simulations. The solid lines are the approximations with function Rrand = a ˜ t + ˜b. The inset shows the dependence of the diffusion coefficient of Pt atoms on the inverse temperature of Cu(001) substrate. The solid line in the inset is the approximation with the functions ˜ B T ). D = D0 exp(−E/k

19

˜ B T ), where E˜ = 1.005 eV and D0 = 2.2 × with function D = D0 exp(−E/k 14 ˚2 10 A /s. We have shown that two processes simultaneously occur after the initial expansion of the Pt cluster in the Cu(001) surface: (i) the linear decrease of the cluster dense core and the loose structure, and (ii) the Brownian motion of Pt atoms on the periphery of the cluster (where alloy has a disordered structure). The following relations were obtained: E2 > E3 > E˜ and E˜ ≈ 1 eV. These results are easy to understand. Indeed, the activation energy for the jump of a Pt atom in the topmost layer of the Cu(001) substrate is the sum of the Gibbs free energy of the formation of a vacancy in the first layer GF = 0.393 eV and the value of the diffusion barrier of this jump (Ed2 = 0.535 eV, see Table 2): EA = GF + Ed2 = 0.928 eV. We see that E˜ ≈ EA . Therefore, atoms on the periphery of the Pt cluster move like single atoms which do not interact with other Pt atoms. The presence of a vacancy in the loose structure is energetically unfavourable and the presence of a vacancy in the dense core is even more unfavourable. Therefore, the ˜ following inequalities are hold: E2 > E3 > E. Finally, we note that a simple one-layer model of the dissolution of a Pt cluster in the Cu(001) substrate considered above allows us to make important conclusions about the nature of the dissolution of Pt clusters in the general case. We expect that the dissolution of the three-dimensional cluster after the initial expansion is also the superposition of two simultaneous processes: the linear decrease of the cluster dense core and the Brownian motion of atoms on the periphery of the cluster. 3.3. Electromigration of small vacancy clusters We used the following approximation to include the electromigration in our SLkMC model. Let us consider the diffusion of a single atom. The following term should be added to the values of the diffusion barriers [56] 1 ∗ Cu,P t = − ZCu,P ∆EEM t eρ (jR) , 2

(10)

∗ where ZCu,P t e is the effective charge of the atom, ρ is the bulk resistivity of copper, j is the average electron current density and R is the vector connecting the initial and the final positions of the atom. The coefficient ∗ ρ∗vac = Zvac eρ for vacancies in copper is −6.6 · 10−7 eV · m/A [57]. Thus, ∗ ∗ in all simulations the effective charge ZCu e = −Zvac e = −ρ∗vac /ρ is used. Unfortunately, the effective charge of Pt in Cu(001) surface is unknown.

20

∗ We perform the simulations with different ratios γ = ZP∗ t /ZCu similarly to Ref. [26]. If the unit Cartesian vectors ex and ey are oriented along the [011] and [0¯11] crystallographic directions, then the equation (10) can be written in the form [27, 28] Cu,P t ∆EEM = αCu,P t (∆nx cos φ + ∆ny sin φ) ,

(11)

where αCu,P t = −ρ∗Cu,P t jr0 /2 is the parameter of electromigration, r0 = 2.556 ˚ A is the first nearest neighbour distance, φ is the angle between ex and j, r0 ∆nx = (Rex ) and r0 ∆ny = (Rey ) are the displacements of the atom along the ex and ey directions. The equation (11) can be generalized for the case of the dimer diffusion X Cu,P t dimer ∆EEM (∆nx,i cos φ + ∆ny,i sin φ) , (12) = αi i=1,2

where αiCu,P t is the parameter of electromigration of the i-th atom in the dimer, r0 ∆nx,i and r0 ∆ny,i are the displacements of the i-th atom. The following designations are used below: αCu = α and αP t = γα, γ ∈ [−1, 1]. Now let us discuss the electromigration of a vacancy clusters in the Pt/Cu(001) surface. A vacancy cluster is formed from the randomly distributed vacancies in the absence of an electric field at the preliminary stage of the simulations. The main stage of the simulations begins after the cluster reaches its equilibrium form. Typical current densities in electromigration experiments are j ∼ 109 ÷ 1012 A/m2 [58–60] and the parameter of electromigration has a value of α ∼ 10−7 ÷ 10−4 eV. At the same time, a typical value of the diffusion barriers for atomic jumps in the topmost layer of the Cu(001) surface is of the order of 0.1 eV. Therefore, it is difficult to simulate the electromigration of vacancy clusters at α < 10−4 eV. Fortunately, the drift velocity of vacancy clusters is directly proportional to the parameter α over a wide range of α up to α = 10−3 eV [24, 27]. Thus, all further calculations will be performed with the parameter of electromigration α = 1 meV. The electromigration of vacancy clusters was simulated with different values of the parameter γ in the range of [−1, 1]. We found that the drift velocity of vacancy clusters does not depend on the parameter γ within the limits of errors. The structure of the surface alloy also does not depend on the parameter γ. Therefore, the averaged values of the drift velocity of vacancy clusters obtained for different values of the parameter γ are presented below. 21

Figure 7: The normalized drift velocity of a vacancy cluster in the topmost layer of the Cu(001) substrate as the function of the concentration of Pt atoms for the vacancy clusters consisting of 10, 40, and 80 vacancies. The substrate temperature is 300 K and α = 1 meV. The sample standard deviations are presented in the plot. The analogous dependence for the Co/Cu(001) alloy from [27] is presented in the inset.

22

The main features of the electromigration of vacancy clusters obtained on the clean Cu(001) surface [28] (linear drift velocity dependence on the parameter α, oscillatory dependence on the cluster size, independence on the angle φ) remain valid also in the case of the Pt/Cu(001) surface alloy. At the same time, the existence of the embedded Pt atoms leads to the decrease of the drift velocity of vacancy clusters. To discuss this effect, it is convenient to define the normalized drift velocity as Vn = V /V0 , where V and V0 are the drift velocities of vacancy clusters in the Pt/Cu(001) alloy and in the clean Cu(001) surface at the same conditions. The dependencies of the normalized drift velocity on the concentration of Pt atoms at 300 K are presented in Fig. 7. The reason for the decrease of the drift velocity of vacancy clusters is the following. The energy increases by 0.148 eV when a vacancy approaches to the Pt atom. Therefore, repulsive forces act between vacancy clusters and Pt atoms. Thus, vacancy clusters avoid getting close to the Pt atoms. Figure 7 shows that the normalized drift velocity of vacancy clusters does not depend on their size. This situation differs significantly from the electromigration of the vacancy clusters in the Cu(001) surface with the embedded Co clusters [27]. In the case of Co clusters the normalized drift velocity of the vacancy clusters dramatically decreases at a certain concentration of the embedded Co atoms nCo (see inset in Fig. 7). This concentration is different for the vacancy clusters of different sizes and can be determined from the condition Dvac (N ) ≈ ∆(nCo ), where Dvac (N ) is the diameter of the vacancy cluster consisting of N atoms, and ∆(nCo ) is the average distance between the Co clusters. The embedded Co clusters are stable and immobile at room temperature [51]. Thus, they play the role of immobile point defects. At the same time, Pt atoms can move in the Cu(001) surface. Therefore, vacancy clusters move in the topmost layer of the Pt/Cu(001) surface almost without deformations. The drift velocity of vacancy clusters decreases mainly because the diffusion rate of Pt atoms is limited. The diffusion rate of Pt atoms almost does not depend on the size of a vacancy cluster. Thus, the normalized drift velocity of vacancy clusters is, indeed, independent on their size. 4. Conclusion In summary, we have investigated the following diffusion-mediated processes in Pt/Cu(001) surface alloy: (i) order-disorder phase transition, (ii) dissolution of small Pt clusters, and (iii) electromigration of small vacancy 23

clusters in the topmost layer. Let us outline the main results of these investigations. Pt atoms embed to the topmost layer of the Cu(001) substrate after their deposition and the Pt/Cu(001) surface alloy is formed. The structure of this alloy depends on the concentration of Pt atoms, the evolution time and the temperature of the substrate. Deposition of small amount of Pt atoms (nP t < 0.1 ML) results in the formation of a disordered alloy. The short-range order, corresponding to the c(2 × 2) structure, was observed at higher concentrations of Pt atoms. The formation of the short-range order was discussed in detail at the concentration of Pt atoms nP t = 0.25 ML. At higher concentration of Pt atoms (nP t ∈ [0.3, 0.5] ML) the long-range order may occur in the alloy. The relaxation times depend exponentially on the inverse temperature and the concentration of Pt atoms. The long-range order was not observed at room temperature. We considered the dissolution of round single-layer clusters in the topmost layer of the Cu(001). The following two processes occur simultaneously after the initial expansion of the cluster dense core: the linear decrease of the cluster core and the Brownian motion of the Pt atoms at the periphery of the cluster. The dissolution times depend exponentially on the initial size of the cluster. Our results were obtained for a simple 2D model but we expect the similar behaviour in the case of the 3D Pt clusters embedded into Cu substrate. Finally, we investigated electromigration of small vacancy clusters. Embedded Pt atoms significantly decrease the mobility of vacancy islands. We emphasize that dependence of the electromigration rate on the impurity concentration essentially depends on whether this impurity forms compact clusters or not. We discussed this issue in detail on the example of the Pt/Cu(001) and Co/Cu(001) systems. We had to fit the universal parameters of the interatomic TB-SMA potentials for the Pt–Cu system in order to perform the SLkMC simulations. Special attention was paid to the quality of reproduction of the binding energies of atoms and the diffusion barriers in the fitting procedure. We hope that the fitted parameters of the TB-SMA potentials will be useful in various investigations related to the diffusion of atoms in the Pt/Cu(001) and Pt/Cu(111) surface alloys.

24

Acknowledgements The research is carried out using the equipment of the shared research facilities of HPC computing resources at Lomonosov Moscow State University [61, 62]. The work was supported by the Foundation for the Advancement of Theoretical Physics and Mathematics “BASIS”. References [1] J. Knudsen, A. U. Nilekar, R. T. Vang, J. Schnadt, E. L. Kunkes, J. A. Dumesic, M. Mavrikakis, F. Besenbacher, A Cu/Pt Near-Surface Alloy for Water-Gas Shift Catalysis, J. Am. Chem. Soc. 129 (2007) 6485–6490. [2] Z. He, Z. Wang, X. Liu, P. Han, Preparation of TiAl–Cr surface alloy by plasma-surface alloying technique, Vacuum 89 (2013) 280 – 284. [3] J. Greeley, M. Mavrikakis, Near-surface alloys for hydrogen fuel cell applications, Catal. Today 111 (2006) 52 – 58. [4] G. W. Graham, P. J. Schmitz, P. A. Thiel, Growth of Rh, Pd, and Pt films on Cu(100), Phys. Rev. B 41 (1990) 3353–3359. [5] R. Belkhou, J. Thiele, C. Guillot, Growth of PtCu(100): formation of a surface alloy, Surf. Sci. 377-379 (1997) 948 – 952. [6] M. Walker, C. Parkinson, M. Draxler, C. McConville, Growth of thin platinum films on Cu(100): CAICISS, XPS and LEED studies, Surf. Sci. 584 (2005) 153 – 160. [7] Y. Shen, J. Yao, D. O’Connor, B. King, R. MacDonald, A search for clock reconstruction in fcc (001) surfaces induced by monolayer metal films: PdCu(001), PtCu(001) and Pd/Pt/Cu(001), Solid State Commun. 100 (1996) 21 – 26. [8] Y. Shen, D. O’Connor, K. Wandelt, Composition and structure of Cu3 Pt(001): a (1x1) Cu termination with c(2x2) underlayer ordering, Surf. Sci. 406 (1998) 23 – 31. [9] E. AlShamaileh, H. Younis, C. Barnes, K. Pussi, M. Lindroos, A tensor LEED determination of the structure and compositional profile of a Cu(100)-c(2x2)-Pt surface alloy, Surf. Sci. 515 (2002) 94 – 102. 25

[10] E. AlShamaileh, C. J. Barnes, A. Wander, Cu-capped surface alloys of Pt/Cu{100}, J. Phys.: Condens. Matter 15 (2003) 1879–1887. [11] J. P. Reilly, D. O’Connell, C. J. Barnes, Modification of formate stability by alloying: the Cu(100)-c(2x2)-Pt system, J. Phys.: Condens. Matter 11 (1999) 8417–8430. [12] G. Demarco, J. E. Garcs, G. Bozzolo, Growth of Pt/Cu(100): an atomistic modeling comparison with the Pd/Cu(100) surface alloy, Surf. Sci. 526 (2003) 309 – 322. [13] G. Pal, G. Lefkidis, W. Hbner, Ab initio Investigation of Pt Dimers on Cu(001) Surface, J. Phys. Chem. A 113 (2009) 12071–12078. [14] G. Pal, G. Lefkidis, W. Hbner, Electronic excitations and optical spectra of Pt2 and Pt4 on Cu(001) modeled by a cluster, Phys. Status Solidi (b) 247 (2010) 1109–1115. [15] I. M. Goldby, B. von Issendorff, L. Kuipers, R. E. Palmer, Gas condensation source for production and deposition of size-selected metal clusters, Rev. Sci. Instrum. 68 (1997) 3327–3334. [16] K. Wegner, W. J. Stark, S. E. Pratsinis, Flame-nozzle synthesis of nanoparticles with closely controlled size, morphology and crystallinity, Mater. Lett. 55 (2002) 318 – 321. [17] C. G. Zimmermann, M. Yeadon, K. Nordlund, J. M. Gibson, R. S. Averback, U. Herr, K. Samwer, Burrowing of Co Nanoparticles on Clean Cu and Ag Surfaces, Phys. Rev. Lett. 83 (1999) 1163–1166. [18] D. De Salvador, E. Napolitani, G. Bisognin, A. Carnera, E. Bruno, S. Mirabella, G. Impellizzeri, F. Priolo, Experimental evidences for two paths in the dissolution process of B clusters in crystalline Si, Appl. Phys. Lett. 87 (2005) 221902. [19] K. Masuo, O. Plaksin, Y. Fudamoto, N. Okubo, Y. Takeda, N. Kishimoto, Effects of laser irradiation on nanoparticle evolution in SiO2 implanted with Cu ions, Nucl. Instrum. Methods Phys. Res., Sect. B 247 (2006) 268 – 270.

26

[20] N. Castin, M. I. Pascuet, L. Malerba, Modeling the first stages of Cu precipitation in α-Fe using a hybrid atomistic kinetic Monte Carlo approach, J. Chem. Phys. 135 (2011) 064502. [21] M. J. Attardo, R. Rosenberg, Electromigration Damage in Aluminum Film Conductors, J. Appl. Phys. 41 (1970) 2381–2386. [22] J. R. Black, Electromigration - A brief survey and some recent results, IEEE Trans. Electron Devices 16 (1969) 338–347. [23] P. S. Ho, T. Kwok, Electromigration in metals, Rep. Prog. Phys. 52 (1989) 301–348. [24] H. Mehl, O. Biham, O. Millo, M. Karimi, Electromigration-induced flow of islands and voids on the Cu(001) surface, Phys. Rev. B 61 (2000) 4975–4982. [25] A. Latz, S. P. Sindermann, L. Brendel, G. Dumpich, F.-J. M. zu Heringdorf, D. E. Wolf, Anisotropy of electromigration-induced void and island drift, J. Phys.: Condens. Matter 26 (2013) 055005. [26] M. Jongmanns, A. Latz, D. E. Wolf, Impurity-induced island pinning during electromigration, Europhys. Lett. 110 (2015) 16001. [27] S. V. Kolesnikov, A. M. Saletsky, Kinetic Monte Carlo simulation of small vacancy clusters electromigration on clean and defective Cu(100) surface, Eur. Phys. J. B 92 (2019) 14. [28] S. V. Kolesnikov, A. M. Saletsky, Electromigration of Small Vacancy Clusters on the (100) Copper Surface, JETP Lett. 108 (2018) 18–22. [29] R. A. Friesner, Ab initio quantum chemistry: Methodology and applications, PNAS 102 (2005) 6648–6653. [30] M. Ben-Nun, T. J. Martnez, Ab Initio Quantum Molecular Dynamics, John Wiley & Sons, Ltd, 2002, pp. 439–512. [31] J. Hafner, Ab-initio simulations of materials using VASP: Densityfunctional theory and beyond, J. Comput. Chem. 29 (2008) 2044–2078. [32] W. Kohn, Nobel Lecture: Electronic structure of matter—wave functions and density functionals, Rev. Mod. Phys. 71 (1999) 1253–1266. 27

[33] F. Cleri, V. Rosato, Tight-Binding Potentials for Transition Metals and Alloys, Phys. Rev. B 48 (1993) 22–33. [34] V. Rosato, M. Guillope, B. Legrand, Thermodynamical and Structural Properties of f.c.c. Transition Metals Using a Simple Tight-Binding Model, Philos. Mag. A 59 (1989) 321–336. [35] N. A. Levanov, V. S. Stepanyuk, W. Hergert, D. I. Bazhanov, P. H. Dederichs, A. Katsnelson, C. Massobrio, Energetics of Co adatoms on the Cu(001) surface, Phys. Rev. B 61 (2000) 2230–2234. [36] N. N. Negulyaev, V. S. Stepanyuk, W. Hergert, P. Bruno, J. Kirschner, Atomic-scale self-organization of Fe nanostripes on stepped Cu(111) surfaces: Molecular dynamics and kinetic Monte Carlo simulations, Phys. Rev. B 77 (2008) 085430. [37] S. A. Dokukin, S. V. Kolesnikov, A. M. Saletsky, A. L. Klavsyuk, Growth of the Pt/Cu(111) surface alloy: Self-learning kinetic Monte Carlo simulations, J. Alloys Compd. 763 (2018) 719–727. [38] A. F. Voter, A Monte Carlo Method for Determining Free-Energy Differences and Transition State Theory Rate Constants, J. Chem. Phys. 82 (1985) 1890–1899. [39] O. Trushin, A. Karim, A. Kara, T. S. Rahman, Self-Learning Kinetic Monte Carlo Method: Application to Cu(111), Phys. Rev. B 72 (2005) 115401. [40] G. Nandipati, N. Govind, A. Andersen, A. Rohatgi, Self-Learning Kinetic Monte Carlo Simulations of Al Diffusion in Mg, J. Phys.: Condens. Matt. 28 (2016) 155001. [41] A. Latz, L. Brendel, D. E. Wolf, A Three-Dimensional Self-Learning Kinetic Monte Carlo Model: Application to Ag(111), J. Phys.: Condens. Matt. 24 (2012) 485005. [42] F. El-Mellouhi, N. Mousseau, L. J. Lewis, Kinetic Activation-Relaxation Technique: An off-Lattice Self-Learning Kinetic Monte Carlo Algorithm, Phys. Rev. B 78 (2008) 153202.

28

[43] S. Kolesnikov, A. Klavsyuk, A. Saletsky, Self-Organisation and Magnetic Properties of Co Nanostructures Embedded in a Cu(100) Surface, Surf. Sci. 612 (2013) 48 – 56. [44] B. R. Brooks, R. E. Bruccoleri, B. D. Olafson, D. J. States, S. Swaminathan, M. Karplus, CHARMM: A Program for Macromolecular Energy, Minimization, and Dynamics Calculations, J. Comput. Chem. 4 (1983) 187–217. [45] N. N. Negulyaev, V. S. Stepanyuk, P. Bruno, L. Diekh¨oner, P. Wahl, K. Kern, Bilayer Growth of Nanoscale Co Islands on Cu(111), Phys. Rev. B 77 (2008) 125437. [46] O. V. Lysenko, V. S. Stepanyuk, W. Hergert, J. Kirschner, Mesoscopic Relaxation in Homoepitaxial Metal Growth, Phys. Rev. Lett. 89 (2002) 126102. [47] S. V. Kolesnikov, A. L. Klavsyuk, A. M. Saletsky, Simulation of the formation of vacancies upon scanning of Cu(100) surface, JETP Lett. 89 (2009) 471–474. [48] T. Siahaan, O. Kurnosikov, H. J. M. Swagten, B. Koopmans, S. V. Kolesnikov, A. M. Saletsky, A. L. Klavsyuk, Co Diffusion in the NearSurface Region of Cu, Phys. Rev. B 94 (2016) 195435. [49] V. S. Stepanyuk, A. L. Klavsyuk, W. Hergert, A. M. Saletsky, P. Bruno, I. Mertig, Magnetism and Structure of Atomic-Size Nanocontacts, Phys. Rev. B 70 (2004) 195420. [50] S. Dokukin, S. Kolesnikov, A. Saletsky, Dendritic growth of the Pt-Cu islands on Cu(111) surface: Self-learning kinetic Monte Carlo simulations, Surf. Sci. 689 (2019) 121464. [51] S. Kolesnikov, A. Klavsyuk, A. Saletsky, The role of the diffusion of dimers in the formation of Co nanostructures embedded into Cu(100) surface, Eur. Phys. J. B 86 (2013) 399. [52] G. Henkelman, H. J´onsson, Improved tangent estimate in the nudged elastic band method for finding minimum energy paths and saddle points, J. Chem. Phys. 113 (2000) 9978–9985. 29

[53] G. Marsaglia, A. Zaman, Some portable very-long-period random number generators, Comput. Phys. 8 (1994) 117–121. [54] J. M. Ziman, Models of disorder : the theoretical physics of homogeneously disordered systems, Cambridge University Press, 1979. [55] L. D. Landau, E. M. Lifshitz, Course of Theoretical Physics: Volume 5, Statistical Physics, 3rd ed., Butterworth-Heinemann, 2013. [56] R. S. Sorbello, Theory of Electromigration, in: Solid State Physics, volume 51 of Solid State Physics, Academic Press, 1998, pp. 159 – 231. [57] J. Hoekstra, A. P. Sutton, T. N. Todorov, A. P. Horsfield, Electromigration of vacancies in copper, Phys. Rev. B 62 (2000) 8568–8571. [58] C. Tao, W. G. Cullen, E. D. Williams, Visualizing the Electron Scattering Force in Nanostructures, Science 328 (2010) 736–740. [59] R. Hoffmann-Vogel, Electromigration and the structure of metallic nanocontacts, Appl. Phys. Rev. 4 (2017) 031302. [60] A. Latz, S. Sindermann, L. Brendel, G. Dumpich, F.-J. Meyer zu Heringdorf, D. E. Wolf, Simulation of electromigration effects on voids in monocrystalline Ag films, Phys. Rev. B 85 (2012) 035449. [61] V. Sadovnichy, A. Tikhonravov, V. Voevodin, V. Opanasenko, ”Lomonosov”: Supercomputing at Moscow State University, in: Contemporary High Performance Computing: From Petascale toward Exascale, Chapman & Hall/CRC Computational Science, Boca Raton, United States, Boca Raton, United States, 2013, pp. 283–307. [62] V. Voevodin, A. Antonov, D. Nikitenko, P. Shvets, S. Sobolev, I. Sidorov, K. Stefanov, V. Voevodin, S. Zhumatiy, Supercomputer lomonosov-2: Large scale, deep monitoring and fine analytics for the user community, Supercomput. Front. Innov. 6 (2019) 4–11.

30

We have no conflict of interest to declare.