University of Birmingham DFT global optimisation of gas-phase and MgO-supported sub-nanometre AuPd clusters

a The Birmingham Parallel Genetic Algorithm (BPGA) has been adopted for the global optimization of free and MgO(100)-supported Pd, Au and AuPd nanocluster structures, over the size range N = 4–10. Structures were evaluated directly using density functional theory, which has allowed the identification of Pd, Au and AuPd global minima. The energetics, structures, and tendency of segregation have been evaluated by different stability criteria such as binding energy, excess energy, second difference in energy, and adsorption energy. The ability of the approach in searching for putative global minimum has been assessed against a systematic homotop search method, which shows a high degree of success.


Introduction
Nanomaterials have at least one dimension on the nanometer scale (1-100 nm). They have recently emerged as new materials that bridge the gap between atoms or molecules and bulk materials, and they have attracted remarkable interest owing to their numerous potential applications. 1 Nanomaterials contribute to many new technological applications in various fields, such as medicine, materials, physics, and chemistry. These new applications came as a consequence of their novel chemical and physical properties which are due to electronic and quantum effects and the high surface-area-to-volume ratio. [1][2][3] The ability to control the surface and structural properties of nanostructures in the nanometre range and their suitable integration with different scientific research concepts have attracted widespread interest from researchers because of their use in many emerging technological applications. [4][5][6][7] AuPd nanostructures, in particular, have been investigated previously for a number of applications, including catalytic applications. AuPd catalysts have been found to be promising candidates for a wide variety of chemical reactions, such as cyclohexane oxidation, 8 NO reduction, 9 CO oxidation, 10,11 direct synthesis of hydrogen peroxide, 12,13 and synthesis of aldehydes from primary alcohols. 14 Theoretically, AuPd systems have been studied to rationalize their catalytic activities. Studies have been conducted on some of the AuPd clusters in the gas phase 15 and supported on surfaces (e.g. MgO(100) and TiO 2 (110) slabs) [16][17][18] in addition to the simulation of the interaction of pure Au and Pd nanoparticles and AuPd nanoparticles with atoms or small molecules such as S, H, NO, and CO. 9,10,19,20 However to date there have been few studies of small AuPd clusters supported on MgO. Subnanometre clusters are groups or aggregates of a few to tens of metal atoms which are (o1.0 nm) in size. These smaller clusters are of interest in catalysis due to the potential for enhanced activity and selectivity.
The rarity of theoretical studies concerned with N = 4-10 AuPd nanoalloys, whether in the gas phase or supported on MgO, has led us to select the Au-Pd system to be studied. In order to find the global minima for the pure and mixed AuPd clusters in the gas phase and on MgO, we have used the Birmingham Parallel Genetic Algorithm (BPGA) to enable global optimisation directly at the DFT level. 21

BPGA-DFT
The BPGA-DFT approach was applied for cluster sizes ranging from N = 4-10 for all compositions of free and MgO(100)supported AuPd nanoalloys, as well as the pure Au and Pd clusters. Gamma-point DFT calculations were performed with the Vienna ab initio Simulation Package (VASP) code 22 utilising projected-augmented wave (PAW) pseudopotentials and the PBE exchange correlation functional. 23,24 A plane-wave basis set was used. The energy was truncated at 400 eV. Methfessel-Paxton smearing, with a sigma value of 0.01 eV, was implemented to improve convergence. 25 The Birmingham Parallel Genetic Algorithm (BPGA) has been adopted for the evaluation of potential cluster structures. This method is the latest open-source genetic algorithm 26 improving on the Birmingham Cluster Genetic Algorithm (BCGA), a genetic algorithm for determining the lowest energy isomers of nanoparticles and nanoalloys up to approximately 100 atoms. 27 BPGA employs a pool methodology to evaluate structures in parallel. In each run, multiple BPGA instances are implemented, and in each instance, a set of processes are run in parallel and independently. 15,28 Numerous random geometries are initially generated to form a population. 21 The generated structures of a given population are then geometrically relaxed (local energy minimized). 15 Once the local minimization of the initial pool structures is completed, crossover and mutation are conducted on the lowest energy individuals in the population.
Clusters are selected for crossover using tournament selection with crossover being performed using the cut-and-splice method introduced by Deaven and Ho. 29 The mutation operators are set as a homotop-swap for nanoalloys and a random atom displacement for the pure clusters. The energy of any newly created structure is compared with those of the other structures in the pool. The highest-energy isomer is replaced with the new lower energy isomer.
For the supported clusters, owing to the high computational cost of MgO(100) slab relaxation, the slab is not relaxed during the local minimization, but the cluster geometries are optimised in the presence of the fixed slab. The supported cluster is optimised within a sphere placed at 1.5 Å over a 6 Â 6 Â 2 slab of MgO(100) with a 14.7 Å vacuum spacing (as shown in Fig. 1), to ensure that there are no cluster-cluster or other non-physical interactions arising due to the periodic boundary conditions. The interactions between the cluster and the surface have been replicated using two layers of MgO(100). 17,[30][31][32] The efficacy of replicating a surface behaviour and the cluster properties by using two layers of the MgO slab has previously been confirmed. 33,34 Random rotation of the cluster regarding the fixed surface is used as a mutation operator when optimising the supported clusters.

Energy calculations
The average binding energy per atom E b is given by where m, n, E Au , and E Pd are the numbers of Au and Pd atoms and the electronic energies of a single Au or Pd atom, respectively, and N is the total number of atoms (N = m + n).
The stability of each cluster, relative to its neighbours, is indicated by the second difference in energy D 2 E which is given by where A is Au or Pd, E (A N ) corresponds to the total energy of the N-atom cluster and E (A (N+1) ) and E (A (NÀ1) ) are the neighbouring clusters, with one atom more and one atom less, respectively. Evaluation of the effect of mixing in binary nanoalloys has been achieved by calculating the excess energy D which is given by where E (Au m Pd n ) is the total energy of the nanoalloy and E (Au N ) and E (Pd N ) are the energies of the pure Au and Pd clusters with the same total number of atoms as Au m Pd n . The adsorption energy, E ads , of the Au m Pd n cluster on the MgO(100) support was calculated by where E (slab+Au m Pd n ) is the energy of the MgO(100)-supported cluster, E (slab) is the energy of the MgO(100) surface, and E (Au m Pd n ) is the energy of the free Au m Pd n cluster, locally minimized in the gas phase.

Results and discussion
3.1 Global optimisation of free and supported Au, Pd, and AuPd clusters 3.1.1 Gold clusters. The putative global minima for pure Au clusters, 4 r N r 10, are shown in Fig. 2 and their energies, structures, and point groups are listed in Table S1 (see the ESI †). Putative global minima for the MgO(100)-supported Au clusters, 4 r N r 10, are shown in Fig. 3    The best structures obtained for all studied free Au clusters in this size range show planar configurations, as previously reported. [35][36][37] Our results satisfactorily concur with those presented by previous theoretical studies on Au clusters up to 10 atoms [38][39][40] and by experimental research on Au cation clusters up to 7 atoms. 41 In addition, our results agree with the findings of Zhao and co-workers 42 concerning the planar structures of small Au clusters up to 6 atoms. However, some differences are observed for Au 7 , Au 8 , and Au 9 , which are predicted to be 3D by Zhao. 42 Global minima for the supported Au clusters are all found to be still planar, with some slight deviation from the planarity for Au 9 and Au 10 . They are found to lie roughly perpendicular to the MgO surface, due to the ''metal-on-top'' effect. 32 Bonding of two atoms of Au to O atoms has led to an increase in the Au-Au distance from (2.68) to (3.76) Å, forming an elongated Au 4 cluster on the surface. 43 N = (5-8) Au clusters on the surface have similar structures to that in the gas phase, whereas Au 9 adopts a different pseudo-planar configuration and Au 10 keeps the gas phase structure but with a bend in one of its edges.
The preference for Au clusters with planar structures may be attributed to the non-additive many-body interactions in comparison with the additive two-body forces in Au atoms, 44 and the involvement of the d electrons of Au in bonding in planar structures is higher than in 3D ones. 39 Some researchers have considered ''relativistic effects'' 45 in interpreting the preference for planar structures, attributing such a preference to the decrease in the 5d-6s orbital spacing that strengthens s-d hybridization. 46 3.1.2 Palladium clusters. Fig. 5 and 6 show the putative global minima for free and supported pure Pd clusters, and Table S1 lists the energies, structures, and point groups for free clusters (see the ESI †).
In contrast to Au clusters, Pd clusters do not adopt 2D structures; instead, they all have 3D motifs that mostly favour    deltahedral (triangular faced) compact structures as lowest-energy configurations. Global minima for the gas phase Pd clusters are all found to be different from their supported structures with the exception of Pd 4 which remains a tetrahedron.
Landman et al. 47 investigated the structural properties of neutral Pd N isomers with N = 1-7 using DFT. These researchers clarified that the global minima are 3D for clusters with more than 3 atoms, which is in full agreement with the gas-phase Pd clusters for sizes 4 r N r 7 reported here. Thus, the global minimum of free Pd 5 is the trigonal bipyramid. On the surface this structure distorts to give a square-based pyramidal structure. The four Pd atoms of the square base are bonded to four different O atoms of the slab. The octahedral structure of free Pd 6 alters to form a bicapped tetrahedron on the surface, forming three Pd-O bonds.
We have previously presented BCGA-DFT studies of free Pd 8 , Pd 9 , and Pd 10 clusters, 16,48 whose putative global minima were identified as a dodecahedron, an icosahedral fragment, and an incomplete centered icosahedron, respectively. These observations conform to our findings for free clusters. Our supported structures of Pd 7 , Pd 8 , and Pd 9 are a capped trigonal prism, a distorted square antiprism, and a distorted tricapped octahedron, respectively. The structure of Pd 10 is an icosahedral fragment in the gas-phase. On the surface this global minimum distorts to give a more complex fused structure.
3.1.3 Gold-palladium clusters. The putative global minima for all compositions of free and supported Au m Pd n clusters, 4 r (m + n) r 10, are shown in Fig. 7-16 and the energies, structures, and point groups of all free clusters are given in Table S2 (see the ESI †). On going from monometallic clusters to nanoalloys, there is an increase in difficulty of identifying the global minimum (GM) structure due to the presence of ''homotops'' in addition to the size effect. 49 All the favoured structures of Au 1 Pd n n = 3-9 compositions are 3D. For N = 4-7 and N = 9, doping one Au atom into the pure Pd clusters yields geometries which are similar to the pure Pd clusters, as previously reported. 35,50 Au 1 Pd n n = 7 and 9 are a bicapped octahedron and a capped face-sharing octahedron, respectively. The Au atom occupies a low-connectivity vertex or an edge position, as previously reported for larger clusters. 51 On the surface, Au 1 Pd n structures are all found to be different from their gas-phase counterparts, with the exception of Au 1 Pd 3 , which remains a tetrahedron.
The 2D structures of Au clusters, discussed in Section (3.1.1), remain the global minima when they are doped with a single Pd atom. The Pd atoms are located in or close to the centre of the clusters. Au m Pd 1 m = 3-5 clusters have the same structures as the corresponding pure Au clusters, as previously reported. 35,50 Au m Pd 1 m = 6 adopts a Pd-centred planar hexagon which is also the core of the clusters with 7 and 8 Au atoms. The cluster for m = 9 also has a planar global minimum. 2D-3D structural transitions take place at sizes N = 7-10 on the surface.
It is clear that Pd-rich nanoalloys tend to adopt 3D geometries while 2D-3D structural transitions occur for gold clusters when adding more than one Pd atom. It is worthwhile mentioning that most of the lowest energy structures of free Au m Pd n clusters for (n + m) = 4-10, which were calculated recently by Zanti and Peeters 35 and Palagin and Doye, 50 are in good agreement with our observations.
The supported and gas-phase structures of Au 1 Pd 3 and Au 2 Pd 2 are both tetrahedra. For Au 3 Pd 1 , the supported structure is similar to that of the free cluster, both having a planar is planar. For Au 1 Pd 4 , the supported structure is similar to that of supported Pd 5 (discussed above), having a squarebased pyramid structure. For Au 2 Pd 3 , the supported cluster is a distorted homotop of the gas phase structure, which allows      all 3 Pd atoms to bind to the surface. A structural transition occurs for the gas-phase Au 3 Pd 2 cluster from 3D to 2D when supported on the surface, having a similar structure to the supported Au 4 Pd 1 structure.
The gas-phase Au 1 Pd 5 structure is an octahedron and Au 2 Pd 4 and Au 3 Pd 3 are both tetrahedra whereas Au 4 Pd 2 and Au 5 Pd 1 are an edge-bridged-capped trigonal bipyramid and planar triangle, respectively. Structural changes occur for all of these global minima on the surface except for Au 5 Pd 1 . The interactions between Au and Pd atoms and the MgO(100) surface favours face-capped square-based pyramid and facecapped trigonal bipyramid structures for supported Au 1 Pd 5 and Au 2 Pd 4 clusters, respectively. The supported configurations of Au 3 Pd 3 , Au 4 Pd 2 and Au 5 Pd 1 are all 2D structures. For Au 3 Pd 3 on the surface, the three Pd atoms form a linear chain bonding to O atoms of the surface, resulting in a planar parallelogram. DFT local minimisation using the VASP code was used to confirm the parallelogram configuration by assessing an alternative planar triangular structure (as previously found for Au 3 Ir 3 43 ) for this composition (see Fig. 12). The minimisation showed that the putative GM parallelogram is indeed lower in energy than the triangle for this composition. However, both Au 4 Pd 2 and Au 5 Pd 1 are found to adopt the planar triangular structure on the surface with clear preference for Pd bonding to O atoms on the surface.
A capped octahedron is the global minimum for both Au 1 Pd 6 and Au 2 Pd 5 , which distort on the surface to form a   bicapped square pyramid and a square pyramid fused with a trigonal bipyramid, respectively, bound to the surface by Pd-O bonds. In the gas phase, the GM for Au 3 Pd 4 , Au 4 Pd 3 and Au 5 Pd 2 all have edge-bridged bicapped tetrahedral structures, with the Pd atoms occupying the central tetrahedral core. On MgO, they all form low-symmetry fused structures. Au 6 Pd 1 has a Pd-centred planar hexagonal structure in the gas phase, while on the surface a complex 3D structure is preferred.
Both Au 1 Pd 7 and Au 2 Pd 6 have the same bicapped octahedral structure (with the caps adopting a ''para'' position) in the gasphase and both adopt structures similar to Au 2 Pd 5 , but with an additional atom capping the square pyramid. The gas-phase Au 3 Pd 5 and Au 4 Pd 4 clusters both have a ''meta-'' bicapped octahedral structure, but on the MgO surface Au 3 Pd 5 adopts a polytetrahedral geometry, while Au 4 Pd 4 keeps the bicapped octahedral structure, but adopts a different homotop so that 3 of the Pd atoms are in contact with the surface. The edgebridged pentagonal bipyramid structure of Au 5 Pd 3 distorts to a polytetrahedral structure on the surface. The gas-phase Au 6 Pd 2 cluster is a hexagonal bipyramid and its global minimum on the surface is similar to that of the supported Au 5 Pd 2 structure with an additional face cap. Homotop swap and structural changes occur for the edge-bridged planar hexagon structure, the global minimum of gas-phase Au 7 Pd 1 , giving a 3-dimensional bent triangle bonded by the single Pd atom and 4 Au atoms to the O atoms of the MgO slab.
The gas-phase Au 1 Pd 8 global minimum is a face-sharing bioctahedral structure. Au 2 Pd 7 -Au 5 Pd 4 have structures that can be described as the result of fusing an octahedron and a trigonal bipyramid, while Au 6 Pd 3 and Au 7 Pd 2 (bicapped pentagonal bipyramid) have structures which are fragments of the centred icosahedron. Au 8 Pd 1 has a planar, bi-edge-bridged centred hexagon. On the MgO surface, complex fused structures are found, mostly consisting of 2-close packed layers of metal atoms, with Pd enrichment closest to the MgO surface.
Doping up to three Au atoms changes the structure of gas-phase Pd 10 (discussed earlier) from an icosahedral fragment to two face-sharing octahedra, with a capping atom at the joint. The tetrahedral tetracapped octahedron is the GM for Au 4 Pd 6 -Au 8 Pd 2 , while Au 9 Pd 1 has a planar, bihexagonal structure. On the MgO(100)-support, bilayer close packed structures predominate, again with the MgO-cluster interface enriched in Pd to maximise Pd-O bonding.

Energetic analysis
To evaluate the stability of free and supported clusters and predict the structural preferences for magic sizes, we have to calculate the excess energy D, the second difference in energy D 2 E, the binding energy E b , and the adsorption energy E ads , which are defined in eqn (1)-(4). Tables S4 and S5 list the values of these energies for free and surface supported clusters (see the ESI †).
Looking at the elemental properties of Au and Pd (shown in Table 1), 52 such as surface energy (E sur ), cohesive energy (E coh ), atomic radius (r a ), and electronegativity (w), we can simply predict the type of segregation in AuPd nanoalloys.
The larger values of cohesive energy and surface energy of Pd compared with Au explain the tendency of Pd atoms to occupy the centres and cores of clusters and the preference of Au atoms for surface positions. The atomic radius of Pd is smaller than for Au. This also favours the aggregation of Pd atoms in the cores of clusters, whereas electron transfer between Pd and Au atoms supports the mixing of Au and Pd, but this is a weak effect as the electronegativities of Au and Pd are quite similar.
The effect of mixing in a cluster system can be studied by calculating the excess energy, D. Excess energy plots for free AuPd clusters (m + n) = 4-10 are shown in Fig. 17 and 18. For (m + n = 4) clusters, Au 1 Pd 3 shows a strong demixing tendency, whereas Au 3 Pd 1 exhibits the best mixing among clusters of this size. All (m + n = 5) clusters favour mixing and the strongest mixing is for Au 2 Pd 3 . However, all mixed (m + n = 6) clusters are energetically unfavourable relative to pure Pd and Au clusters. For (m + n = 7) clusters, all clusters show favourable mixing, with the largest tendency for Au 6 Pd 1 . For (m + n = 8) clusters, only Au 1 Pd 7 shows unfavourable mixing and the most favourable mixing is for Au 3 Pd 5 . The clusters of (m + n = 9 and 10) favour mixing, with the most favourable mixing being for Au 8 Pd 1 and Au 6 Pd 4 .
The mixing energy values for clusters with the same size show that variations in composition play a more important role than geometric effects in determining cluster stabilities.
The relative stabilities of clusters can be studied by calculating the second difference in energy D 2 E, which indicates the stability of an N-atom cluster with respect to neighbouring sizes. Fig. 19 shows a plot of the second difference in energy D 2 E for free and supported Au and Pd clusters. The most relatively stable clusters are indicated by significant positive peaks in D 2 E.
The results show that free Au 6 (2D) and Pd 6 (3D) clusters have a high relative stability compared to their neighbours. Accordingly, N = 6 represents a magic size for both 2D Au and 3D Pd structures. The magic size N = 6 for both free Au and Pd clusters explains the positive values of mixing energy of nanoalloys compared to monometallic clusters at this size. On the MgO surface, Au 6 is still found to have the highest relative stability due to it retaining the same triangular structure both when free and supported on MgO. However, Pd clusters show a noticeable (stable) positive peak in D 2 E for the supported cluster with N = 5, rather than 6. This is because the structure of Pd 6 changes from an octahedron to a bicapped tetrahedron on the MgO surface, in addition to the cluster-substrate interaction effects mentioned in Section 3.1.2.
The relative stabilities of nanoclusters can be obtained by calculating the binding energy per atom, E b . A plot of the binding energies for studied clusters is shown in Fig. 20. The binding energy increases with increasing cluster size, as previously reported. 35 The plot also illustrates the relative stability of nanoalloys compared to the same size of free Au and Pd clusters, with one exception for N = 6, which reinforces the ''magic size'' hypothesis mentioned above.
There are differences between the energies of the Mg(100)supported global minima and the energies of the supported clusters with the slab removed. Tables S4 and S5 (see the ESI †) list the adsorption energies E ads of the supported clusters and Fig. 21 shows a graphical representation. The plots reveal that the Pd and Pd-rich clusters tend to have the most negative E ads as there are a greater number of strong Pd-O interactions.
Comparing the E ads values for the clusters as a function of size is more complicated, however, due to the interplay between    structural changes and the number of interfacial Pd atoms; as for smaller clusters a higher proportion of atoms are in contact with the oxide substrate. The higher adsorption energies of Pd-doped Au clusters can contribute to improve catalyst performance by decreasing cluster diffusion rates and, hence, suppressing sintering.

Systematic homotop search
The global minima (GM) of AuPd nanoalloys have many symmetry inequivalent homotops which may have been missed by the BPGA-DFT search. Based on the GM of monosubstituted clusters and some supported-clusters, the structural energy for each symmetry inequivalent homotop was studied by DFT local minimisation using the VASP code. Tables S6-S9 (see the ESI †) show energies and structures of homotops for free monosubstituted clusters and some MgO(100)-supported clusters. Fig. 22 and 23 show plots of the relative energy DE against symmetry inequivalent homotop structures for these clusters.
The BPGA-DFT search for all free monosubstituted clusters N = 4-10, successfully found the lowest energy homotop as the global minimum. For N = 4, Au 1 Pd 3 does not have any symmetry inequivalent homotops whereas Au 3 Pd 1 has one symmetry inequivalent homotop with a relative energy of 0.50 eV.
For N = 5, the only symmetry inequivalent homotop for Au 1 Pd 4 is 1.04 eV higher in energy than the GM found by BPGA-DFT. Au 4 Pd 1 has two symmetry inequivalent homotops and both are higher in energy than the GM, by 0.92 and 0.96 eV, as shown in Fig. 22(a).
For N = 6, Au 1 Pd 5 does not have any symmetry inequivalent homotops whereas Au 5 Pd 1 has just one symmetry inequivalent homotop which is 2.86 eV less stable than the GM. For N = 7, Au 1 Pd 6 has two symmetry inequivalent homotops with relative energies of 0.11 eV and 1.93 eV, as shown in Fig. 22(b), whereas Au 6 Pd 1 has one symmetry inequivalent homotop with a relative energy of 0.28 eV. For N = 8, Au 1 Pd 7 has just one symmetry inequivalent homotop with a relative energy of 0.35 eV, whereas Au 7 Pd 1 has four symmetry inequivalent homotops which can clearly be seen in Fig. 22(c). For N = 9, Au 1 Pd 8 also has one symmetry inequivalent homotop with a relative energy of 0.15 eV, whereas Au 8 Pd 1 has three symmetry inequivalent homotops which can be seen in Fig. 22(d). Finally, Fig. 22(e) and (f) shows the relative energies for symmetry inequivalent homotops of Au 1 Pd 9 and Au 9 Pd 1 .
The BPGA-DFT search for MgO(100)-supported clusters generally successfully finds the energy homotop as the global minimum. Fig. 23, for example, shows the relative energy for symmetry inequivalent homotop structures of the supported Au 5 Pd 1 cluster.

Conclusions
The use of the BPGA-DFT approach has successfully allowed the global optimization of free and MgO(100)-supported N = 4-10 Au, Pd, and AuPd clusters. Significant structural differences between the gas-phase and surface-supported are revealed by global optimization in the presence of a MgO slab in addition to the clear homotop swap behaviour in mixed AuPd clusters, in order to increase the number of stabilising Pd-O bonds. For singly doped clusters, the GM structures were confirmed by homotop reminimisation.
The BPGA-DFT calculations show that the Pd and Pd-rich clusters prefer 3D structures, while structural transitions occur for Au configurations from 2D to 3D upon adding more than one Pd atom. The planar structures of free Au global minima are all found to remain planar on the MgO(100) surface, with  some deviation from the planarity for supported Au 9 and Au 10 . They are also found to lie roughly perpendicular to the MgO surface, due to the ''metal-on-top'' effect. Global minima for the gas phase Pd clusters are all found to be different from their supported structures with the exception of Pd 4 , which retains its tetrahedral structure.
Mixing energies show the strong tendency of free Au-Pd clusters to alloy, with the exception of N = 6, which confirms the magic size hypothesis of pure Pd and Au clusters at this size (confirmed by findings of the binding energy E b for nanoalloys compared to pure clusters and the second difference in energy D 2 E calculations for pure Au and Pd clusters).
The adsorption energy E ads results reveal that the Pd and Pd-rich clusters tend to have the most negative E ads values. This could improve catalytic performance by suppressing cluster sintering.