Stability and mobility of supported Ni_n (n = 1–10) clusters on ZrO₂(111) and YSZ(111) surfaces: a density functional theory study

Abstract

The performance of supported metal catalysts, such as nickel nanoparticles decorating yttria-stabilized zirconia (YSZ), depends on their microstructure and the metal–support interface. Here, we have used spin polarized density functional theory (DFT) to evaluate different Ni cluster geometries and determined the electronic structure of the most stable configurations. We have described the interaction of Ni_n (n = 1–10) clusters supported on the cubic ZrO₂(111) and YSZ(111) surfaces, which show a preference for pyramidal shapes rather than flat structures wetting the surface. The interfacial interaction is characterized by charge transfer from the cluster to the surface. We also show how yttrium, present in YSZ, affects the Ni–Ni interaction. Through analysing the difference between the cohesive energy and the clustering energy, we show the preference of Ni–Ni bond formation over Ni-surface interaction; this energy difference decreases with the increase of the Ni_n cluster size. From the evaluation of the Ni atomic hopping rates on YSZ, we have demonstrated that under different temperature conditions, Ni atoms aggregate with other atoms and clusters, which affects the cluster size stability.

---

1. Introduction

There is a growing interest in the stabilities and properties of supported clusters and metal-oxide interfaces, since many industrial processes use oxide-supported metal particles as catalysts, i.e. in the water gas shift reaction,¹ steam reforming processes,^2,3 microelectronics, sensors and solid oxide fuel cells (SOFC).^4–6 These processes depend, among other things, on the particle sizes and the particle–support interactions.^7–9 For instance, the performance of the electrode in SOFCs – electrochemical devices that convert chemical energy into electrical energy – depends on the microstructure and the size distribution of the metal particles and the oxide phase in the cermet¹⁰ (ceramic matrix containing metal nanoparticles). A common electrode in SOFCs is metallic nickel, used as an electron conductor, supported on an oxygen conductor, usually an oxide. For instance, Ni supported on zirconia (ZrO₂) doped with yttria (Y₂O₃) (Ni/YSZ), is a suitable material for SOFCs as it has a high mixed electronic-ionic conductivity and high mechanical strength.¹¹ The key reactions in these devices occur at the triple phase boundary (TPB) where the gas phase, Ni particles and the YSZ surface meet. Therefore, tuning and stabilizing the Ni microstructure will affect the SOFC performance and electrode lifetime.¹² For example, the performance decay of the SOFC is associated with the sintering of Ni particles through atomic or cluster diffusion across the oxide support.¹³

While atomic-level information on the metal-oxide interface is difficult to obtain experimentally, computational tools facilitate the evaluation of individual factors and thereby can provide crucial insight.^14,15 For example, previous density functional theory (DFT) calculations have been used to study the interaction of certain metals on γ-Al₂O₃;^16–19 Rajesh et al. studied the adsorption of Au_n (n = 2–10) on α-Al₂O₃ (0001), showing favourable metal clustering rather than wetting i.e. spreading across the surface.²⁰ Di Valentin et al. investigated the stability and growth of Ni clusters on the MgO surface,²¹ where small Ni_n (n ≤ 9) clusters adsorb weakly on the surface and their diffusion can be stopped by surface defects such as oxygen vacancies. Molina et al.²² have shown that an Au₂₀ cluster keeps its tetrahedral geometry upon adsorption on the MgO(100) surface. Here, we have systematically analysed the interaction of Ni_n (n = 2–10) clusters on cubic zirconia (c-ZrO₂)(111) and YSZ(111) surfaces in order to characterize the electronic structure of the interface and the cluster geometry as a function of their size, and to evaluate the stability and sintering rate of Ni clusters on both c-ZrO₂(111) and YSZ(111) surfaces.

2. Models and computational methods

We have employed calculations based on density functional theory (DFT)²³ as implemented in the Vienna Ab initio Simulation Package (VASP),^24–27 where the generalized gradient approximation (GGA)²⁸ with the Perdew–Burke–Ernzerhof (PBE) density functional was used to approximate the exchange–correlation functional. Long-range dispersion interactions were described by using the semi-empirical method of Grimme²⁹ and we have employed the projected augmented wave method (PAW)³⁰ to describe the interaction between the valence and the core electrons, where O (2s, 2p), Ni (3d, 4s), Zr (4d, 5s) and Y (4s, 4p, 4d, 5s) atomic orbitals have been treated as valence electrons. All the calculations were spin polarized. The unit cell of c-ZrO₂ has the fluorite structure with space group Fm3m. We have converged the bulk energy for different k-points (from 7 × 7 × 7 to 11 × 11 × 11) and cutoff energies (from 350 eV to 550 eV) with a convergence criterion of 0.001 eV per atom. The bulk convergence was reached for a kinetic energy cutoff of 500 eV and a Monkhorst–Pack grid of 9 × 9 × 9 mesh of k-points. For the 1 × 1 slab calculations, we have tried four sets of k-points (from 6 × 6 × 1 to 9 × 9 × 1) with the same convergence criterion and convergence was reached for a 7 × 7 × 1 mesh of k-points. For the 2 × 2 slab calculations, three sets of k-points were tried (3 × 3 × 1, 4 × 4 × 1, and 5 × 5 × 1) with convergence reached for the 3 × 3 × 1 mesh. The tetrahedron method with Blöchl corrections was used to define the Brillouin zone.

The ZrO₂(111) surface was obtained with the METADISE code,³¹ which takes into account the atomic charges and the periodicity of the plane leading to a stacking of atomic layers with zero dipole moment perpendicular to the surface plane.

We have modelled the most stable zirconia surface, ZrO₂(111)^32–35 (Fig. 1(a) and (b)), considering a slab with nine periodic atomic layers (three O–Zr–O trilayers). The five top atomic layers were relaxed during the optimization, while the four bottom layers were kept fixed at their bulk equilibrium position to represent the rest of the crystal. The vacuum size between slabs was set to 15 Å in order to minimise interactions with the perpendicular images.

Fig. 1 (a) and (b) represent the side and top view, respectively, of the oxygen terminated ZrO₂(111) surface. (c) and (d) correspond to the side and top view, respectively, of the YSZ(111) surface. O_u corresponds to the topmost oxygen atom, while O_d is located in the lower oxygen layer of the uppermost ZrO₂ trilayer. Color key: O, Zr and Y are respectively represented by red, grey and cyan spheres. The yellow spheres correspond to oxygen vacancies.

In our previous investigation,³⁵ we determined the band gap to be 3.14 for the ZrO₂(111) surface. Although the calculated band gap is in agreement with previous computational studies,^36,37 it is lower than the experimentally measured bulk values (5–6 eV),^38,39 as DFT calculations tend to underestimate the band gap.

We substituted two Zr with two Y atoms and removed one oxygen atom in order to obtain the yttria-stabilized zirconia (YSZ) surface with a dopant concentration of 9.09%, which is in the range of an optimized oxygen transport YSZ electrode (8–10%). We have previously evaluated the energy of different non-equivalent configurations of the YSZ(111) surface and found that the most stable structure had the two yttrium atoms located in the top and sub O–Zr–O trilayers.³⁵ The oxygen vacancy is positioned in the lower oxygen layer of the top O–Zr–O trilayer and at the next nearest neighbour site (NNN) to the two Y atoms (Fig. 1(c) and (d)).

There are many possible shapes and adsorption sites for Ni_n clusters. We have chosen the Ni_n geometries based on the stability of fcc Ni surfaces^40–43 and small particles⁴⁴ to form triangular based structures in agreement with previous studies. We have therefore grown clusters exposing the Ni(111)-facet^45,46 on the ZrO₂(111) and YSZ(111) surfaces, where the initial atomic spacing between the Ni atoms is 2.2 Å, similar to the interatomic distance in the Ni(111) surface. We have based the initial positions of our clusters on our previous findings.³⁵ On the ZrO₂(111) surface, the clusters interact with the oxygen atoms located in the top oxygen layer, O_u (Fig. 1(a) and (b)), which act as a trap for the Ni. On the YSZ(111) surface, Ni atoms adsorb preferentially on top of the oxygen vacancy, away from the Y atoms,³⁵ and hence, we have grown the Ni_n (n = 2–10) clusters away from the Y atoms.

We have evaluated the stability of the growing Ni cluster by using the clustering energy per atom (E_clus), eqn (1):

where E_Nin/surf, E_surf and E_Ni are the energies of the slab with the Ni_n cluster, the clean oxide surface and the Ni metal atom in the bulk, respectively, and n is the number of adsorbed Ni atoms. A positive clustering energy means that the Ni cluster growth is favourable. The evaluation of the cohesive energy (E_coh) (eqn (2)) allowed us to understand the interactions between the Ni atoms in the clusters:where E_Nin is the energy of the Ni_n cluster with the same adsorbed geometry in vacuum (single point calculation). Finally, we have studied the perpendicular interaction between the surface and the Ni clusters by determining the interaction energy (E_int) (eqn (3)):where is the energy of the geometry of the clean oxide surface taken from the optimized Ni_n/surf structure and n_C is the number of Ni atoms in contact with the surface. The charge transfer between the surface and the molecules was analysed using Bader analysis as implemented in the Henkelman algorithm.⁴⁷

We also implemented the DFT results in kinetic Monte Carlo simulations to study the coalescence of Ni atoms and clusters on YSZ(111), thus providing information regarding coverage and time at specific temperatures.

3. Results and discussion

3.1. Structural analysis

3.1.1. Ni_n/ZrO₂(111). In our previous investigation,⁴⁸ we showed that Ni nanoclusters adopt a hexagonal arrangement similar to the shape of the (111)-facet of face-centred cubic metals. The clusters sit on the surfaces, optimizing the interactions with O_u atoms and transferring some electrons to the surface. Based on these preliminary findings, we have positioned up to Ni₁₀ by maximizing the interactions with O_u surface atoms, where we have considered flat and pyramidal morphologies. The Ni_n (n = 1–5) atoms prefer a 3D structure rather than wetting the surface with planar sheets.⁴⁸ For example, we have tried four initial configurations for the Ni₆ cluster on ZrO₂(111) (Fig. S1 in the ESI†): two flat and two in a three-dimensional shape. The most stable shape has five Ni atoms at the base of the structure, see Fig. 2, maximising their interaction with the O_u atoms and having an E_clus of 1.51 eV. The average Ni–Ni distance in this cluster is 2.4 Å, which is 0.2 Å shorter than those in the Ni(111) surface. The optimized structure is different from the initial guess, where before relaxation the Ni atoms were positioned in a flat distribution in contact with the surface. After relaxation, the atoms adopted a Ni(111)-facetted tetrahedron shape on the zirconia surface. Similar studies^20,49 have shown the preference of aggregation of metal clusters on oxide surfaces. For example, Zhang et al.⁴⁹ have shown that for n ≥ 4 Ru_n adopts a 3-D geometry on TiO₂(101), whereas Rajesh et al.²⁰ have shown that the aggregation of Au atoms is more favourable than wetting the Al₂O₃ surface.

Fig. 2 Side views of the most stable configurations of Ni_n/ZrO₂(111) systems. Colour key: red, grey and blue spheres correspond to oxygen, Zr and Ni atoms, respectively.

The analysis of the interfacial distances in the Ni₆/ZrO₂(111) system shows that each Ni atom is at ∼1.9 Å from its nearest O_u neighbour. In addition, we note that these O_u atoms are shifted from their initial position: the O_u–O_u distance is 3.6 Å in the clean ZrO₂(111) surface, whereas it is approximately 4.1 Å upon Ni₆ adsorption. The O_u atoms are pushed apart in order to accommodate the cluster and minimise the forces on it. The adsorption of this cluster also affects the electronic structure of ZrO₂. We have evaluated the atomic charges (Table 1) and observed a total charge transfer of 0.5 e from the metal cluster to the surface. The charge transfer is slightly higher for the Ni₉ and Ni₁₀ clusters: +0.6 and +0.7 e, respectively.

Table 1 Calculated clustering (E_clus), cohesive (E_coh) and perpendicular interaction (E_int) energies (in eV), and Ni_n charges (e) of the most stable configuration of both Ni_n/ZrO₂(111) and Ni_n/YSZ(111) systems. Ni_n/ZrO₂(111) (n = 2–5) values are taken from ref. 48 for comparison

Nin	E clus		E coh		E int		Total Nin charge
	ZrO2(111)	YSZ(111)	ZrO2(111)	YSZ(111)	ZrO2(111)	YSZ(111)	ZrO2(111)	YSZ(111)
Ni2	2.16	2.39	3.83	3.85	−2.38	−1.64	+0.2	+0.1
Ni3	1.77	2.02	3.38	3.38	−2.16	−1.53	+0.3	+0.2
Ni4	1.57	1.80	3.11	3.06	−2.04	−1.82	+0.4	+0.1
Ni5	1.53	1.79	2.88	2.89	−2.05	−1.60	+0.4	+0.2
Ni6	1.51	1.76	2.77	2.67	−1.94	−1.51	+0.5	+0.2
Ni7	1.47	1.67	2.67	2.54	−1.95	−1.42	+0.5	+0.2
Pyr-Ni8	1.42	1.68	2.51	2.51	−1.95	−1.50	+0.5	+0.2
Flat-Ni8	—	1.67	—	2.66	—	−1.66	—	+0.3
Ni9	1.39	1.57	2.40	2.32	−2.01	−1.31	+0.6	+0.3
Ni10	1.31	1.58	2.30	2.22	−1.87	−1.26	+0.7	+0.2

---

This charge transfer causes an electronic rearrangement of the surface atoms interacting with the cluster: the negative charge of the O_u atoms decreases from −1.2 to −1.1 e, while the charge of the Zr atoms (near the cluster) decreases from 2.3 to 2.1 e. We have also noted that only the five Ni atoms directly bonded to the surface transfer charge to ZrO₂(111), since each of those Ni atoms are +0.1 e while the Ni atom at the top of the cluster remains uncharged. Generally, Ni atoms at the vertex are either negatively charged (−0.1 e) or neutral, which indicates that the atoms at the base of the cluster transfer charge to the surface, while the low coordinated metal atoms, i.e. vertex and corners, accumulate electron density.⁴⁴ Thus, the Ni located at the vertex of the cluster could be a source of electrons for an eventual reaction with electron receptor molecules approaching the cluster from the gas phase.⁵⁰ Analysing the electron density difference between the cluster and the surface (Fig. 3), we have confirmed a charge accumulation between the Ni atoms directly bonded to the surface and the surface atoms; there is no charge relocation between the top Ni atom which remains fully metallic. We also note the well-localised accumulation of electron density between the O_u and the Ni atoms, indicating the formation of a Ni–O bond. A previous report has also shown this electron rearrangement between Pt clusters and the zirconia surface.⁵¹

Fig. 3 Calculated electron density difference between Ni_n clusters and the surface for the most stable configuration of the Ni_n/ZrO₂(111) systems. Colour key: red, grey and blue spheres correspond to oxygen, Zr and Ni atoms, respectively. Isosurface level: 0.005 electrons per ų.

Next, we have analysed the interaction of Ni₇, Ni₈, Ni₉ and Ni₁₀ with the ZrO₂(111) surface where we have tried, respectively, five, six, three and four non-equivalent initial configurations (Fig. S2 to S5 in the ESI†). The most stable shapes are shown in Fig. 2 and the calculated clustering energies are 1.47, 1.42, 1.39 and 1.31 eV, respectively (Table 1). The four configurations have a similar pyramid shape and the only difference is the number of the Ni atoms in the top Ni layer of each cluster (Fig. 2).

3.1.2. Ni_n/YSZ(111). In the YSZ system, we have again built several clusters with non-equivalent initial adsorption configurations, both three-dimensional and flat, and adsorbed them on top of the YSZ(111) surface. A single Ni atom on top of the YSZ(111) surface sits on top of the oxygen vacancy, as far away as possible from the yttrium atom, with a clustering energy of 2.35 eV.³⁵

We have optimised several initial configurations for the flat Ni₂/YSZ(111) and Ni₃/YSZ(111) clusters (Fig. S6 and S7 in the ESI†) and the most stable ones are shown in Fig. 4, with an average clustering energy of 2.21 eV (Table 1).

Fig. 4 Side views of the most stable configurations of Ni_n/YSZ(111) systems. Colour key: red, grey, blue and cyan spheres correspond to oxygen, Zr, Ni and Y atoms, respectively. The oxygen vacancy is represented by a yellow sphere.

For Ni_n (n = 4–7), the average E_clus is 1.76 eV and the clusters adopt a Ni(111) facet shape, in agreement with a previous study by Li et al.¹⁴ who showed that the Ni₄ cluster prefers to adopt a similar pyramid shape on top of γ-Al₂O₃. Similarly, Carrasco et al.⁵² have also found that the pyramid Ni₄ geometry is preferred over the planar one, when the cluster is adsorbed on the isostructural CeO₂(111) surface.

In the Ni_n/YSZ(111) (n = 2–7) configurations, the Ni-surface atomic spacing is 2.3 Å, i.e. slightly larger than in the pure ZrO₂(111) systems. In addition, each O_u is pushed towards the neighbouring vacancy, since we observe an average decrease of 0.1 Å of the O_u-vacancy distance. This movement explains the preference of the Ni clusters to adsorb on this site: it is easier to drive the O_u towards a vacancy where it tends to fill this defect. Moreover, the most stable adsorption site is the one involving two O_u atoms with one neighbouring oxygen vacancy.

Charge analysis (Table 1) shows a slight charge transfer from the cluster to the surface of ∼+0.2 e. As was observed for the Ni_n/ZrO₂(111) systems, this charge transfer affects the electronic structure of the surface: the charges of Zr and Y atoms (near the metal cluster) decrease from 2.2 to 2.1 e and the negative charge of O_u decreases from −1.2 to −1.1 e. We also observe charge accumulation on the Ni atom at the vertex of the pyramid of the Ni_n/YSZ(111) (n = 4–7) configurations, since this latter atom has a −0.1 e charge. The electron density difference plot (Fig. 5) shows this charge accumulation between the cluster and the surface and shows the localized orbitals of O_u and Y atoms interacting with the cluster.

Fig. 5 Calculated electron density difference for the most stable configuration of the Ni_n/YSZ(111) systems. Colour key: red, grey, blue and cyan spheres correspond to oxygen, Zr, Ni and Y atoms, respectively. The oxygen vacancy is represented by a yellow sphere. Isosurface level: 0.002 electrons per ų.

The flat shapes of the Ni_n/YSZ(111) (n = 2–7) structures, shown in the ESI,† are less stable by ∼0.11 eV (which corresponds to a thermal energy of 638 K), i.e. the Ni atoms prefer to aggregate rather than spread over the surface, as was also observed for the Ni_n/ZrO₂(111) systems.

It is worth noting that for the Ni₈/YSZ(111) cluster, the flat and pyramid shapes are in equilibrium since their clustering energies differ by only 0.01 eV (58 K) (Table 1). These two configurations are also the lowest energy ones from the six configurations considered (Fig. S16 and S17 in the ESI†). In the flat configuration (Fig. S16 in the ESI†) the Ni₈ cluster interacts with six O_u atoms at an average Ni–O_u distance of 2.0 Å. In addition, the oxygen vacancy is filled by one displaced O_u which contributes to the stability of the system. Indeed, the movement of O_u towards the vacancy allows the cluster to optimise its interaction with the other six O_u surface atoms. This could be a drawback, for example in SOFCs, since the role of the oxygen vacancies is to enhance the oxygen transport. As for the three-dimensional shape (Fig. 4), the cluster interacts with five O_u atoms but none of the oxygen vacancies are filled by an O_u atom. This shape is similar to the one found for Ni₇/YSZ(111), although the next Ni atom avoids interaction with Y from the surface and forms a (111)-facetted shape with the rest of the cluster. Bader charge analysis (Table 1) shows that the charge transfer from the cluster to the surface is slightly higher in the flat configuration (+0.3 e) than in the pyramid shape (+0.2 e) owing to the former’s interaction with an extra O_u atom. This comparison shows the importance of the two parameters responsible for the stability of the cluster on top of the surface: the number of O_u atoms interacting with the cluster and the shape of the cluster.

The last cluster studied in this work is Ni₁₀/YSZ(111), where three configurations have been considered (Fig. S20 and S21 in the ESI†), and the most stable shape found is the pyramid one (Fig. 5). The calculated clustering energy is 1.58 eV (Table 1) which is similar to the value found for Ni₉/YSZ(111). This is predictable since both Ni₉/YSZ(111) and Ni₁₀/YSZ(111) have exactly the same shape and the only difference is the tenth Ni atom located at the apex of the pyramid in Ni₁₀/YSZ(111). Therefore, the modification of the surface geometry upon adsorption and the number of O_u atoms involved in the interaction are similar to Ni₉/YSZ(111). As to the Bader charge (Table 1), the cluster is +0.3 e charged and we note that the Ni atom located at the apex of the cluster is −0.2 e charged, while all the other Ni atoms are positively charged (average of +0.1 e). We also note the same variation of the charge of the surface atoms interacting with Ni₁₀ as was observed for Ni₉/YSZ(111). Thus, the Ni atoms close to the surface transfer charge to the surface atoms and to the Ni atom located at the top of the cluster, which can be seen from the electron density difference plot (Fig. 5), where we have a gain of charge between the cluster and the surface.

In general, in Ni_n/YSZ(111) (n = 2–10), the Ni_n clusters interact with O_u atoms with an average distance of 1.9 Å. The metal clusters transfer charge to the surface, depending on the size, ranging from +0.1 to +0.3 e. The three-dimensional cluster shape is more favourable due to the repulsive interaction between Y and Ni.

3.2. Ni sintering on top of ZrO₂(111) and YSZ(111)

We have shown that Ni_n clusters on both ZrO₂(111) and YSZ(111) surfaces prefer to adopt a three-dimensional structure rather than flat shapes for clusters containing at least 4 atoms. The same conclusion was drawn for a similar system, CeO₂-supported Au nanoparticles,⁵³ where it was shown that planar Au₁₃ on top of the CeO₂ surface is unstable compared to three-dimensional Au₁₃. Pan et al.⁵⁴ have also shown that a Ni₄ cluster prefers to adopt a 3-D pyramid shape on top of the γ-Al₂O₃(110) surface, with a large clustering energy. Giordano et al.⁵⁵ have also demonstrated that the Ni₄ cluster prefers to adopt a tetrahedron shape on top of the MgO(001) surface.

The total energies of four individual Ni atoms (4Ni) compared to the Ni₄ cluster adsorbed on the surface are 3.38 eV and 2.20 eV less stable on the ZrO₂(111) and YSZ(111) surfaces, respectively, thus showing that aggregation of the Ni atoms is clearly preferred energetically over dispersion. We have also compared the total energy of two Ni₄ clusters, separated by approximately 6.0 Å, with the Ni₈/surface system and here we also found that (Ni₄ + Ni₄) is less stable than Ni₈, now by 1.67 eV and 1.11 eV on ZrO₂(111) and YSZ(111), respectively.

Furthermore, in Fig. 6(a) we have plotted the clustering energy (E_clus) as a function of the cluster size, for both ZrO₂(111) and YSZ(111) surfaces. This graph shows a probable aggregation of Ni on both surfaces owing to a thermodynamic driving force as the cluster size increases. The trend in the clustering energy shows that for the same cluster size, the clustering energy is lower on ZrO₂(111), in agreement with the two aggregation examples calculated for Ni₄ and Ni₈ clusters. The Y atoms affect the geometry of the surface and anion rearrangement, making the interaction between the clusters and the O_u surface atoms less favourable, thus enhancing the preference for the formation of Ni clusters.

Fig. 6 (a), (b) and (c) represent, respectively, the plot of the clustering energy (E_clus), the difference between the clustering energy and the cohesive energy (E_coh − E_clus) and the interaction energy (E_int) as a function of the cluster size on top of both ZrO₂(111) and YSZ(111) surfaces.

We have also evaluated the difference between the cohesive energy (E_coh) and the clustering energy (E_clus) as a function of the cluster size. The difference between those energies expresses the trend to form a Ni–Ni bond instead of a Ni-surface bond. Fig. 6(b) shows that for both surfaces, this energy difference decreases with increasing cluster size due to the preference of the Ni–Ni interaction over Ni-surface interactions. The E_coh − E_clus graph also shows that for the same cluster size, the energy difference is lower for YSZ(111) than for ZrO₂(111), indicating that the degree of interaction between Ni atoms is greater for Ni_n clusters on top of YSZ(111).

The interaction energy (E_int) calculated as a function of the cluster size (Fig. 6(c)) confirms the affinity of the Ni_n clusters for the ZrO₂(111) surface over YSZ(111). Note that the Ni cluster interaction energy is more favourable for zirconia than yttria-stabilized zirconia, whose preference is even more striking for larger clusters. This demonstrates that cluster aggregation is more favourable on the YSZ(111) surface than on ZrO₂(111), in good agreement with, for instance, the difference in energy between two Ni₄ clusters and a Ni₈ cluster, which is larger on ZrO₂(111) than YSZ(111).

From the energy profiles in Fig. 7, we note that the Ni₅/YSZ(111) configuration is more stable than (Ni + Ni₄)/YSZ(111), indicating that Ni atoms prefer to aggregate to form larger clusters rather than wetting the surface. This is in good agreement with the graph in Fig. 6, where we have shown a decrease of the E_coh − E_clus difference, indicating the energetic favourability of aggregation of Ni atoms over their dispersion on both surfaces. From Fig. 7, we also observe that the activation energy, i.e. the energy difference between the (Ni + Ni₄)/YSZ(111) structure and the transition state, is ΔE_Ni5/YSZ = 0.46 eV, which is lower than the same Ni collection adsorbed on ZrO₂(111) (ΔE_Ni5/ZrO2 = 0.72 eV).⁴⁸ This difference in activation energies for the addition of a Ni atom to a larger cluster implies that single Ni atoms can more easily join a bigger cluster when those metal atoms are adsorbed on the YSZ(111) surface rather than ZrO₂(111). It indicates that the aggregation of Ni atoms to form larger clusters is facilitated when the zirconia surface is doped with yttria, which is in good agreement with our E_clus and E_coh − E_clus plots in Fig. 6.

Fig. 7 Energy profile showing two transition states: from state (Ni + Ni₄)/YSZ(111) to Ni₅/YSZ(111) and from state (Ni + Ni₁₀)/YSZ(111) to Ni₁₀/YSZ(111). Colour key: red, grey, blue and cyan spheres correspond to oxygen, Zr, Ni and Y atoms, respectively. The oxygen vacancy is represented by a yellow sphere.

From the calculated activation energies, we have evaluated the hopping rate of one Ni atom from state (Ni + Ni_n)/YSZ(111) to state Ni_(n+1)/YSZ(111): k_A⃑B = ν exp(−ΔE/k_BT) where A is either state (Ni + Ni₄)/YSZ(111) or (Ni + Ni₁₀)/YSZ(111) and B is either state Ni₅/YSZ(111) or Ni₁₁/YSZ(111) (Fig. 8). The Boltzmann constant is k_B = 8.6 × 10⁻⁵ eV K⁻¹ and the vibrational frequency, ν, is accepted as 10¹² s⁻¹. We have therefore calculated k_A→B for a range of temperatures corresponding to the working temperature of a SOFC (T = 500–900 °C). Fig. 8 shows the variation of k_A→B (for five Ni atoms on the surface) from 1.05 × 10⁸ to 1.09 × 10¹⁰ s⁻¹, which is higher than the values found for ZrO₂(111) in our previous investigations (k_A→B varies from 1.87 × 10⁷ to 7.66 × 10⁸ s⁻¹),⁴⁸ indicating that the diffusion of a Ni atom towards a cluster is more favourable on the YSZ(111) surface. This hopping rate is even higher when the cluster diffuses towards a larger cluster, since for eleven Ni atoms k_A→B varies from 586.99 × 10⁷ to 703.93 × 10⁸ s⁻¹ (Fig. 8). Thus, this evaluation of the hopping rate further strengthens the conclusions drawn from Fig. 6: the aggregation of Ni atoms is preferred over dispersion on the surface and this aggregation is more favourable on the YSZ(111) surface.

Fig. 8 Hopping rate for a Ni atom on YSZ(111).

3.3. Microkinetics for Ni_n/YSZ(111)

We have used kinetic Monte Carlo simulations in order to understand how time (t) and Ni coverage (θ) influence the coalescence of Ni clusters on top of the YSZ(111) surface. To achieve our kinetic Monte Carlo study, we have used the same equations as the ones described in the ESI of the investigation by Tafreshi et al.⁵⁶ Additionally, quantum-mechanical tunnelling⁵⁷ has been taken into account. Coverage has been defined in eqn (4):where θ_*(t) is the coverage of free sites and θ_i(t) is the coverage of the Ni_n clusters. The reactions considered are: Ni_n + Ni₁ ↔ Ni_(n+1). We considered as a free site a 2 × 2 clean YSZ(111) surface, while the occupied site is a 2 × 2 YSZ(111) surface with Ni_n clusters adsorbed on top of it.

The first system considered is a YSZ(111) surface with a coverage of 50% of Ni₁ (θ_Ni1 = 50%) as a starting point. In Fig. 9 we show a schematic representation of the initial conditions where one Ni atom covers 50% of the system, i.e., at t = 0.0 seconds, half of the system is considered as a clean 2 × 2 YSZ(111) surface and the other half is Ni₁ supported on YSZ(111).

Fig. 9 Schematic representation of the coverage of the free sites (θ_*) and occupied sites (θ_i) by Ni atoms on the YSZ(111) surface. Here, θ_i = θ_* = 50%.

The coverage (θ_i(t)) of the Ni_n clusters as a function of time for a fixed temperature of 500 K is represented in Fig. 10(a).

Fig. 10 Plot of the coverage, θ_i, of the Ni_n clusters as a function of time. θ_i has been calculated from t = 0.0 s to t = 3600 s but here we plot θ_i up to t = 100 s, as from t = 50 s, the values of θ_i are unchanged. For Ni₈, both pyramid and flat shapes have been considered. (a) and (b) show θ_Ni1 for an initial coverage of 50% and 100%, respectively. (c) and (d) show θ_Ni10 for an initial coverage of 5% and 10%, respectively.

During the first 50 seconds the concentration of Ni₁ drops to 0.0%. In the meantime, θ_Ni4 increases to reach a plateau of 12%. Thus, at T = 500 K, if 50% of the YSZ(111) surface is covered by single Ni atoms, the metal atoms will aggregate spontaneously to generate Ni₄ clusters. This quick coalescence process is in good agreement with the previous section: for high temperatures, the Ni atoms tend to form clusters on top of YSZ(111) rather than wet the surface.

In Fig. 10(b) we plot the coverage as a function of time, for T = 500 K with an initial θ_Ni1 of 100%. The result is similar to the previous initial condition (θ_Ni1 = 50%): after 50 seconds, all of the single Ni atoms aggregate to form Ni₄ pyramids. Indeed, from t = 50 s, θ_Ni4 = 24% and θ_* = 76%. Furthermore, whatever the initial coverage of single Ni atoms is, the aggregation speed is the same as for both initial coverages (θ_Ni1 = 50% and θ_Ni1 = 100%), within 50 seconds Ni₄ pyramid clusters are generated. We have also calculated the reaction rates and noted that the one corresponding to Ni₄ + Ni → Ni₅ has the lowest value. Thus, Ni atoms will aggregate quickly to generate Ni₄ but the formation of larger clusters is negligible.

In Fig. 10(c) we describe the evolution of a Ni₁₀ cluster on YSZ(111) with a starting coverage of 5%. Within 50 seconds Ni₁₀ coverage decreases (3.7%) in order to generate Ni₉ (0.65%) and Ni₁₁ (0.65%). The formation of the Ni₉ clusters is the consequence of the generation of Ni₁₁: one Ni atom, belonging to Ni₁₀, detaches to aggregate with a neighbouring Ni₁₀ cluster giving Ni₉ and Ni₁₁.

In Fig. 10(d), we show θ_Ni10 with an initial coverage of 10%, and the evolution of the cluster size is similar to an initial θ_Ni10 of 5% because the interaction between metal moieties is neglected in the model.

Finally, this kinetic Monte Carlo study shows that Ni_n always tend to form larger clusters; single Ni atoms and Ni₁₀ clusters generate Ni₄ and Ni₁₁, respectively. However, under experimental conditions no specific cluster size will be dominant, as we saw, for instance, that Ni₁₀ can generate Ni₉ and Ni₁₁.

4. Conclusions

We have used spin polarized DFT and kinetic Monte Carlo calculations to study the interaction of Ni_n (n = 1–10) clusters with both ZrO₂(111) and YSZ(111) surfaces. The general trend observed for Ni_n clusters (n = 2–10) adsorbed on the ZrO₂(111) and YSZ(111) surfaces shows the importance of the interaction of Ni with the O_u atoms, as well as the shapes of the adsorbed Ni_n clusters. We have seen that the O_u atoms facilitate the adsorption of the clusters and that these atoms are shifted from their initial position upon adsorption. In Ni_n/YSZ(111) systems, the clusters pushed the O_u atoms towards filling neighbouring vacancies which could be a drawback in, for instance, SOFCs since these vacancies play an important role in oxygen transport. Bader charge analysis of the clusters revealed that there is charge transfer from the cluster to the surface, in particular from the metal atoms at the interface. The distribution of the charge within the cluster is similar for all clusters with a pyramid shape (Ni_n, n = 4–10): the Ni atoms bound to the surface are positively charged and those at the top of the pyramid are either charge-neutral or even some negative charge is accumulated at the apex. In some instances, the Ni atoms located at the top of the clusters have a non-negligible amount of charge, which could play a role in the adsorption of molecules: for example, an electrophilic molecule would adsorb on top of the Ni cluster, rather than at the meeting point between the cluster and the surface. Finally, from calculation of the clustering and cohesive energies and evaluation of the diffusion barriers and hopping rates, we conclude that, on both ZrO₂(111) and YSZ(111) surfaces, the aggregation of the Ni atoms takes place spontaneously, especially on the YSZ(111) surface. Kinetic Monte Carlo simulations showed that Ni atoms aggregate on the YSZ(111) surface once they are exposed to medium temperature: single Ni atoms tend to form Ni₄ clusters, while Ni₁₀ clusters generate Ni₉ and Ni₁₁ clusters. We also noted that the sintering speed does not depend on the initial coverage.

Acknowledgements

We acknowledge the Engineering and Physical Sciences Research Council (Grant no. EP/K001329 and EP/K016288) for funding. Via our membership of the UK’s HPC Materials Chemistry Consortium, which is funded by EPSRC (EP/L000202), this work made use of the facilities of ARCHER, the UK’s national high-performance computing service, which is funded by the Office of Science and Technology through EPSRC’s High End Computing Programme. The authors also acknowledge the use of the UCL@Legion High Performance Computing Facility, and associated support services, in the completion of this work. Finally, the authors acknowledge the use of the IRIDIS High Performance Computing Facility, and associated support services at the University of Southampton, in the completion of this work. Information about the data underpinning the results presented here, including how to access them, can be found in the Cardiff University data catalogue at http://doi.org/10.17035/d.2018.0050609499.

References

1. C. Rhodes, G. J. Hutchings and A. M. Ward, Water-Gas Shift Reaction: Finding the Mecanistic Boundary, Catal. Today, 1995, 23, 43–58.
2. A. Haryanto, S. Fernando, N. Murali and S. Adhikari, Current Status of Hydrogen Production Techniques by Steam Reforming of Ethanol: A Review, Energy Fuels, 2005, 19, 2098–2106.
3. D. J. Moon, Hydrogen Production by Catalytic Reforming of Gaseous Hydrocarbons (Methane & LPG), Catal. Surv. Asia, 2008, 12, 188–202.
4. R. Grau-crespo, N. C. Hernandez, J. F. Sanz and N. H. de Leeuw, Theoretical Investigation of the Deposition of Cu, Ag, and Au Atoms on the ZrO₂(111) Surface, J. Phys. Chem. C, 2007, 2, 10448–10454.
5. C. R. A. Catlow, S. A. French, A. A. Sokol, M. Alfredsson and S. T. Bromley, Understanding the Interface between Oxides and Metals, Faraday Discuss., 2003, 124, 185.
6. W. Wang, Z. Liang, X. Han, J. Chen, C. Xue and H. Zhao, Mechanical and Thermodynamic Properties of ZrO₂ under High-Pressure Phase Transition: A First-Principles Study, J. Alloys Compd., 2015, 622, 504–512.
7. C. R. Henry, Surface Studies of Supported Model Catalysts, Surf. Sci. Rep., 1998, 31, 231–233.
8. D. W. Goodman, Model Studies in Catalysis Using Surface Science Probes, Chem. Rev., 1995, 95, 523–536.
9. J. Y. Park, L. R. Baker and G. A. Somorjai, Role of Hot Electrons and Metal-Oxide Interfaces in Surface Chemistry and Catalytic Reactions, Chem. Rev., 2015, 115, 2781–2817.
10. S. Kim, H. Moon, S. Hyun, J. Moon, J. Kim and H. Lee, Performance and Durability of Ni-Coated YSZ Anodes for Intermediate Temperature Solid Oxide Fuel Cells, Solid State Ionics, 2006, 177, 931–938.
11. E. Drożdż, J. Wyrwa and M. Rękas, Influence of Sintering Temperature and Aging on Properties of Cermet Ni/8YSZ Materials Obtained by Citric Method, Ionics, 2013, 19, 1733–1743.
12. Y. Chen, S. Chen, G. Hackett, H. Finklea, X. Song and K. Gerdes, Microstructural and Chemical Evolution near Anode Triple Phase Boundary in Ni/YSZ Solid Oxide Fuel Cells, Solid State Ionics, 2011, 204–205, 87–90.
13. S. Gao, J. Li and Z. Lin, Theoretical Model for Surface Diffusion Driven Ni-Particle Agglomeration in Anode of Solid Oxide Fuel Cell., J. Power Sources, 2014, 255, 144–150.
14. J. Li, E. Croiset and L. Ricardez-sandoval, Effect of Metal–Support Interface During CH₄ and H₂ Dissociation on Ni/γ-Al₂O₃: A Density Functional Theory Study, J. Phys. Chem. C, 2013, 117, 16907–16920.
15. T. Ogura, K. Nakao, T. Ishimoto and M. Koyama, Computational Study on Impurities Poisoning and Degradation of an SOFC Anode Based on Density Functional Theory, ECS Trans., 2011, 35, 853–858.
16. D. Cao, Y. Li, J. Wang and H. Jiao, Mechanism of γ-Al₂O₃ Support in CO₂ Reforming of CH₄-A Density Functional Theory Study, J. Phys. Chem. C, 2011, 115, 225–233.
17. L. G. V. Briquet, C. R. A. Catlow and S. A. French, Platinum Group Metal Adsorption on Clean and Hydroxylated Corundum Surfaces, J. Phys. Chem. C, 2009, 113, 16747–16756.
18. D. Mei, J. H. Kwak, J. Szanyi, Q. Ge and C. H. F. Peden, Catalyst Size and Morphological Effects on the Interaction of NO₂ with BaO/γ-Al₂O₃ Materials, Catal. Today, 2010, 151, 304–313.
19. Z. Liu, L. Ma and A. S. M. Junaid, NO and NO₂ Adsorption on Al₂O₃ and Ga Modified Al₂O₃ Surfaces: A Density Functional Theory Study, J. Phys. Chem. C, 2010, 114, 4445–4450.
20. C. Rajesh, S. Nigam and C. Majumder, The Structural and Electronic Properties of Au_n Clusters on the α-Al₂O₃(0001) Surface: A First Principles Study, Phys. Chem. Chem. Phys., 2014, 16, 26561–26569.
21. C. Di Valentin, L. Giordano, G. Pacchioni and N. Rösch, Nucleation and Growth of Ni Clusters on Regular Sites and F Centers on the MgO(001) Surface, Surf. Sci., 2003, 522, 175–184.
22. L. Molina and B. Hammer, The Activity of the Tetrahedral Au₂₀ Cluster: Charging and Impurity Effects, J. Catal., 2005, 233, 399–404.
23. P. Hohenberg and W. Kohn, Inhomogeneous Electron Gas, Phys. Rev., 1964, 136, B864–B871.
24. G. Kresse and J. Furthmüller, Efficient Iterative Schemes for Ab Initio Total-Energy Calculations Using a Plane-Wave Basis Set, Phys. Rev. B: Condens. Matter Mater. Phys., 1996, 54, 11169–11186.
25. G. Kresse and J. Furthmüller, Efficiency of Ab-Initio Total Energy Calculations for Metals and Semiconductors Using a Plane-Wave Basis Set, Comput. Mater. Sci., 1996, 6, 15–50.
26. G. Kresse and J. Hafner, Ab Initio Molecular Dynamics for Open-Shell Transition Metals, Phys. Rev. B: Condens. Matter Mater. Phys., 1993, 48, 13115–13118.
27. G. Kresse and J. Hafner, Norm-Conserving and Ultrasoft Pseudopotentials for First-Row and Transition Elements, J. Phys.: Condens. Matter, 1994, 6, 8245–8257.
28. J. P. Perdew, K. Burke and M. Ernzerhof, Generalized Gradient Approximation Made Simple, Phys. Rev. Lett., 1996, 77, 3865–3868.
29. S. Grimme, Semiempirical GGA-Type Density Functional Constructed with a Long-Range Dispersion Correction, J. Comput. Chem., 2006, 27, 1787–1799.
30. P. E. Blöchl, Projector Augmented-Wave Method, Phys. Rev. B, 1994, 50, 17953–17979.
31. G. W. Watson, E. T. Kelsey, N. H. de Leeuw, D. J. Harris and S. C. Parker, Atomistic Simulation of Dislocations, Surfaces and Interfaces in MgO, J. Chem. Soc., Faraday Trans., 1996, 92, 433.
32. A. Christensen and E. Carter, First-Principles Study of the Surfaces of Zirconia, Phys. Rev. B, 1998, 58, 8050–8064.
33. S. Gennard, F. Cora and C. R. A. Catlow, Comparison of the Bulk and Surface Properties of Ceria and Zirconia by Ab Initio Investigations, J. Phys. Chem. B, 1999, 103, 10158–10170.
34. X. Xia, R. Oldman and R. Catlow, Computational Modeling Study of Bulk and Surface of Yttria-Stabilized Cubic Zirconia, Chem. Mater., 2009, 21, 3576–3585.
35. A. Cadi-Essadek, A. Roldan and N. H. de Leeuw, Ni Deposition on Yttria-Stabilized ZrO₂(111) Surfaces: A Density Functional Theory Study, J. Phys. Chem. C, 2015, 119, 6581–6591.
36. M. Shishkin and T. Ziegler, The Oxidation of H₂ and CH₄ on an Oxygen-Enriched Yttria-Stabilized Zirconia Surface: A Theoretical Study Based on Density Functional Theory, J. Phys. Chem. C, 2008, 112, 19662–19669.
37. M. Alfredsson and C. R. A. Catlow, Modelling of Pd and Pt Supported on the {111} and {011} Surfaces of Cubic-ZrO₂, Phys. Chem. Chem. Phys., 2001, 3, 4129–4140.
38. R. H. French, S. J. Glass, F. S. Ohuchi, Y.-N. Xu and W. Y. Ching, Experimental and Theoretical Determination of the Electronic Structure and Optical Properties of Three Phases of ZrO₂, Phys. Rev. B, 1994, 49, 5133–5142.
39. S. Sayan, R. A. Bartynski, X. Zhao, E. P. Gusev, D. Vanderbilt, M. Croft, M. Banaszak Holl and E. Garfunkel, Valence and Conduction Band Offsets of a ZrO₂/SiO_xN_y/n-Si CMOS Gate Stack: A Combined Photoemission and Inverse Photoemission Study, Phys. Status Solidi, 2004, 241, 2246–2252.
40. Y. Liu, Y. Wang and G. Chen, Nucleation and Mobility Model of Ag_n Clusters Adsorbed on Perfect and Oxygen Vacancy MgO Surfaces, J. Mol. Model., 2011, 17, 1061–1068.
41. C. Jung, H. Tsuboi, M. Koyama, M. Kubo, E. Broclawik and A. Miyamoto, Different Support Effect of M/ZrO₂ and M/CeO₂ (M = Pd and Pt) Catalysts on CO Adsorption: A Periodic Density Functional Study, Catal. Today, 2006, 111, 322–327.
42. Z. Li, C. V. Ciobanu, J. Hu, J.-P. Palomares-Báez, J.-L. Rodríguez-López and R. Richards, Experimental and DFT Studies of Gold Nanoparticles Supported on MgO(111) Nano-Sheets and Their Catalytic Activity, Phys. Chem. Chem. Phys., 2011, 13, 2582–2589.
43. X. Ma, Y. Dai, M. Guo, Y. Zhu and B. Huang, Insights into the Adsorption and Energy Transfer of Ag Clusters on the AgCl(100) Surface, Phys. Chem. Chem. Phys., 2013, 15, 8722–8731.
44. A. Roldan, J. J. M. Ricart and F. Illas, Growth and Properties of Au Nanowires, Mol. Simul., 2009, 35, 1051–1056.
45. J. Li, E. Croiset and L. Ricardez-Sandoval, Carbon Clusters on the Ni(111) Surface: A Density Functional Theory Study, Phys. Chem. Chem. Phys., 2014, 16, 2954–2961.
46. S. Hong, Y.-H. Shin and J. Ihm, Crystal Shape of a Nickel Particle Related to Carbon Nanotube Growth, Jpn. J. Appl. Phys., 2002, 41, 6142–6144.
47. R. F. W. A. Bader, Quantum Theory of Molecular Structure and Its Applications, Chem. Rev., 1991, 91, 893–928.
48. A. Cadi-Essadek, A. Roldan and N. H. de Leeuw, Density Functional Theory Study of Ni Clusters Supported on the ZrO₂(111) Surface, Fuel Cells, 2017, 17, 125–131.
49. S. Zhang, C. Li, H. Yan, M. Wei, D. G. Evans and X. Duan, Density Functional Theory Study on the Metal–Support Interaction between Ru Cluster and Anatase TiO₂(101) Surface, J. Phys. Chem. C, 2014, 118, 3514–3522.
50. S. S. Tafreshi, A. Roldan and N. H. de Leeuw, Density Functional Theory Calculations of the Hydrazine Decomposition Mechanism on the Planar and Stepped Cu(111) Surfaces, Phys. Chem. Chem. Phys., 2015, 17, 21533–21546.
51. A. Eichler, Chemical Characterization of a Zirconia-Supported Pt Cluster, Phys. Rev. B, 2005, 71, 125418-1.
52. J. Carrasco, L. Barrio, P. Liu, J. A. Rodriguez and M. Vero, Theoretical Studies of the Adsorption of CO and C on Ni(111) and Ni/CeO₂(111): Evidence of a Strong Metal–Support Interaction, J. Phys. Chem. C, 2013, 2, 8241–8250.
53. H. Y. Kim, H. M. Lee and G. Henkelman, CO Oxidation Mechanism on CeO₂-Supported Au Nanoparticles, J. Am. Chem. Soc., 2012, 134, 1560–1570.
54. Y. Pan, C. Liu and Q. Ge, Effect of Surface Hydroxyls on Selective CO₂ Hydrogenation over Ni₄/γ-Al₂O₃: A Density Functional Theory Study, J. Catal., 2010, 272, 227–234.
55. L. Giordano, G. Pacchioni, A. M. Ferrari, F. Illas and N. Rosch, Electronic Structure and Magnetic Moments of Co₄ and Ni₄ Clusters Supported on the MgO(001) Surface, Surf. Sci., 2001, 473, 213–226.
56. S. S. Tafreshi, A. Roldan and N. H. de Leeuw, Micro-Kinetic Simulations of the Catalytic Decomposition of Hydrazine on the Cu(111) Surface, Faraday Discuss., 2017, 197, 41–57.
57. G. Henkelman, A. Arnaldsson and H. Jónsson, Theoretical Calculations of CH₄ and H₂ Associative Desorption from Ni(111): Could Subsurface Hydrogen Play an Important Role?, J. Chem. Phys., 2006, 124, 1–9.

---

Footnote

† Electronic supplementary information (ESI) available: Less stable configurations and their respective binding energies of Ni clusters on ZrO₂(111) (Fig. S1 to S5) and YSZ(111) (Fig. S6 to S21) surfaces. See DOI: 10.1039/c7fd00217c

---