Structural transformations in Cu, Ag, and Au metal nanoclusters

Metal nanoclusters constitute an important branch of nanotechnology because they exhibit size- and shape-dependent properties. Typically, metal nanoclusters adopt¹ either the non-crystalline icosahedron (Ih) and decahedron (Dh) motifs or the crystalline octahedron (fcc) motif, with the non-crystalline structures being dominant at smaller sizes but becoming unfavorable at large sizes due to the stress contribution to the energy that is proportional to the volume.^2–4 Since properties of technological interest (catalytic, optical, etc.) depend on the cluster structure, it is crucial to understand the equilibrium structures of metal nanoclusters. For this purpose, computer simulations can be very useful. Most of the studies available in the literature focus on finding the global energy minimum at a given size.^3,5–17 Although this information is important, it is limited in the sense that global minima refers to the structures at 0 K. However, metal nanoclusters are expected to be produced and observed at finite temperatures. In addition, various structural motifs coexist^18,19 at a specific size and temperature. Hence, a method to reliably calculate the equilibrium distribution of various structural motifs in the entire temperature range is essential.

One possible approach is the Harmonic Superposition Approximation (HSA),^20,21 which has been used to study Lennard Jones,^22–25 metal,^18,26 and alloy nanoclusters.^27,28 Briefly, in this method, a large number (>10³) of low-lying minima are sampled from the potential energy surface (PES) to construct an approximation of the partition function. Subsequently, the temperature-dependent probability of an isomer is calculated based on the partition function. HSA captures the structural distribution at low temperatures fairly accurately. However, at higher temperatures, HSA becomes progressively erroneous. This stems mainly from the failure to accommodate the anharmonic effects, which become significant at larger temperatures. Another issue is the difficulty in capturing the melting region. In order to reconstruct the melting region, it is necessary to sample the high energy region of the PES, which would require one to collect a prohibitively large number of minima. Due to these constraints, melting cannot be reliably captured using HSA.

Alternatively, to sample the phase space effectively, one can simulate several replicas²⁹ of the system that are at different temperatures and are allowed to exchange configurations at specific intervals according to a Metropolis-like criterion. This method is referred to as replica exchange or parallel tempering. At higher temperatures, the barriers between various structures are easily overcome, ensuring a good sampling at these temperatures. On the other hand, the exchange of configurations allows the high temperature configurations to cascade to lower temperatures and, in the process, to improve the phase space exploration at lower temperatures as well. Both Monte Carlo^30–33 and molecular dynamics^34,35 can be carried out in conjunction with parallel tempering. In parallel tempering Monte Carlo (PTMC), generally, random displacement moves are employed to sample configurations, which reduces the likelihood of inter-motif transition with increasing cluster size.³⁶ In addition, collective atomic rearrangements³⁷ are involved during inter-motif transitions involving metallic clusters, which might not be straightforward to incorporate into Monte Carlo sampling. As a result, in this work, we carry out parallel tempering with molecular dynamics.

Recently, we have proposed a method¹⁹ that combines HSA and parallel tempering leveraging the advantages offered by these two methods to capture the structural distribution in the entire temperature range (0 K to melting). First, we carry out parallel tempering molecular dynamics (PTMD) with several replicas at temperatures ranging from room temperature to beyond melting. A large collection of local minima is sampled during the PTMD simulations, which are then fed into the HSA calculations. This combined method offers several advantages where HSA and PTMD act in a complementary fashion. The conventional HSA calculations require the collection of a large number of local minima, which are obtained using structure optimization methods.^18,27 In our case, the minima are directly obtained from PTMD simulations without the need to explicitly search for them. HSA can capture the low temperature solid–solid transitions, which might prove elusive for PTMD. On the other hand, PTMD captures the high temperatures and the melting regions accurately, where HSA calculations fail. As a result, the low- and high-temperature regions are accurately captured by HSA and PTMD, respectively. In the intermediate temperature range, HSA and PTMD have a good agreement.

In this work, we apply this method to study the size- and system-dependent structural changes with temperature in Cu, Ag, and Au metal nanoclusters. This is crucial information given their strong influence on the properties of metal nanoclusters. For example, the catalytic activity of metal nanoclusters depends on the structure type and size^38–41 due to the wide variety of catalytic sites.⁴² In addition, the catalytic activity can be enhanced by an ensemble of different geometrical structures in comparison to homogeneously shaped structures.⁴³ Hence, it is essential to gather knowledge on the equilibrium structural distribution where various geometrical motifs coexist.

Several theoretical works have calculated the global minimum structures of Cu, Ag, and Au nanoclusters. Grigoryan et al.⁹ calculated the global minima of Cu clusters up to 150 atoms using the embedded atom method (EAM)⁴⁴ and up to 60 atoms using Gupta⁴⁵ and Sutton-Chen⁴⁶ potentials. Highly stable structures occur at the sizes 13, 19, 55, 92, and 147, with all of them having high symmetry icosahedral structures except 92, which is a chiral structure having T point group symmetry. Most of the structures are icosahedra, with the sizes 4, 17, 26, 28, 29, 91–95 having tetrahedral geometry and 75, 78, 81, 101–103 being decahedra. In the case of Ag nanoclusters of sizes larger than 60 atoms, the decahedron is found to be the dominant motif.^7,10,11,15 There are a few exceptions where the truncated octahedron (fcc) and icosahedron (Ih) are the global minima. Due to their strong relativistic effects,⁴⁷ Au nanoclusters exhibit peculiar structures. At sizes smaller than 40 atoms, Au nanoclusters adopt either planar or hollow cage-like geometries.^48–54 In comparison to Cu and Ag, Au disfavors icosahedral structures. At the magic sizes of 55, 147, and 309, the icosahedron is not the global minimum.^14,18,37 This is also evident over a larger size range (up to 1000 atoms).³ However, when the icosahedral structures are observed in Au nanoclusters, for example, at higher temperatures,^55,56 they typically have “rosette”^57,58 defects on the surface. A “rosette” defect appears when a vertex atom is pushed out to form a six-atom ring with the five neighboring surface atoms, leaving behind a vacancy at the vertex position.

Cu, Ag, and Au clusters have also been studied using density functional theory (DFT) calculations. Generally, ideal structures are considered since global minimum search becomes prohibitive at the DFT level for clusters larger than ∼ 50 atoms.⁵⁹ Roldán et al.⁶⁰ carried out structural analysis of several “magic” sized octahedral Cu, Ag, and Au clusters in the range 38–225 atoms and identified a correlation to estimate cohesive energies in a large size range. Similarly, Kiss et al.⁶¹ studied octahedral and icosahedral Ag clusters (consisting of 6–600 atoms) and observed that the cohesive energy is linear with the inverse of cluster size. Oliveira et al.⁶² showed that Ag icosahedra are energetically stable compared to cuboctahedra through density functional tight binding (DFTB) calculations of “magic” clusters in the range 55–561 atoms.

The picture arising from experiments is more complex since, in experiments, it is often difficult to disentangle kinetic effects from equilibrium ones.¹ Electron microscopy has been used to study the structure of metal nanoclusters of varying sizes and temperatures. Langlois et al.⁶³ prepared Cu nanoparticles in a broad size range of 1–12 nm using thermal evaporation. They observed a significant overlap between icosahedra and decahedra at sizes less than 8 nm, beyond which fcc structures were observed. Volk et al.⁶⁴ analyzed Ag clusters with a size <7 nm grown in superfluid He droplets. The smallest particles were fcc, with decahedra at intermediate sizes and icosahedra at large sizes. However, theoretical predictions^2,3 show that icosahedra are energetically favored at smaller sizes, while fcc are favored at larger sizes, and large icosahedra are likely to be due to kinetically trapped growth on top of smaller decahedra.^65,66

Recently, the structural distribution of size-selected Ag clusters centered around 309 atoms was measured,⁶⁷ finding an abundance of fcc structures with very few icosahedra (2%). This is in contrast to the prediction that icosahedra are the dominant motif around the size 309.³ Wells et al.⁶⁸ calculated the proportion of various motifs in Au₅₆₁, Au₇₄₂, and Au₉₂₃. At these sizes, fcc and decahedra make up 70% of the structures, while icosahedra contribute less than 5%. The finite-temperature distribution of Au₅₆₁ was calculated by Foster et al.⁶⁹ in the temperature range 20–500 °C. Again, icosahedra were almost non-existent beyond 100 °C with less than 3%. At temperatures greater than 125 °C, there is an increase in the proportion of decahedra at the expense of fcc structures. The experiments establish a lack of preference for the icosahedral motif in Au nanoclusters, in agreement with the findings of Gupta potential and DFT calculations.⁷⁰ 

From a theoretical and experimental viewpoint, it is essential to have knowledge of the equilibrium proportion of various structural motifs as a function of temperature. In this work, we calculate the structural distribution of Cu, Ag, and Au metal nanoclusters at the sizes 90, 147, and 201, which fall in the size range of 1–2 nm. These were chosen to highlight the size- and system-dependent structural changes. 147 and 201 are “magic” sizes, corresponding to the perfect icosahedron (147) and the regular truncated octahedron (201). It is generally assumed that “magic” sized structures have energetic stability. Our results show that this assumption is not always true. Finally, we chose 90 to look at non-magic sized structures.

We use the tight binding model within the second moment approximation (TBSMA),⁷¹ which is also referred to as the Gupta⁴⁵ potential or Rosato–Guillope–Legrand (RGL)⁷² potential, to model the atom–atom interactions in Cu, Ag, and Au nanoclusters. The parameters of the Gupta potential have been taken from Ref. 2. The interaction potential of Au gives an accurate description of the experimental cluster structures in the gas phase⁶⁸ and on MgO substrates.⁷³ In addition, this potential agrees well with DFT calculations in the prediction of surface “rosette” defects in icosahedra⁵⁷ and the tendency to disfavor icosahedra.⁷⁰ Coming to Ag and Cu, the Gupta potentials correctly predict the stability of Mackay stacking over anti-Mackay stacking in icosahedral clusters in line with the DFT calculations (see supplementary material in Ref. 74). In Ag₅₈₆, the fcc structure is energetically preferred in comparison to the icosahedron, which is also the case according to DFT.⁷⁵ 

Gupta potential predicts correctly that Ag icosahedra are energetically stable compared to cuboctahedra, which agrees with the DFTB calculations⁶² (see the plot of the energy difference between cuboctahedron and icosahedron in supplementary material, Fig. S1). At size 147, icosahedra and decahedra are the prominent motifs. In order to assess the competition between these motifs, we have carried out DFT calculations for Cu₁₄₇ and Ag₁₄₇. For Au₁₄₇ clusters, we refer to the calculations performed previously.⁷⁰ DFT calculations were carried out using the Quantum ESPRESSO⁷⁶ code. Projected augmented wave (PAW)⁷⁷ pseudopotentials were used with the Perdew–Burke–Ernzerhof (PBE)⁷⁸ exchange-correlation functional. An energy cutoff of 45 Ry was used for both Ag and Cu, while a charge density cutoff of 181 Ry and 236 Ry were used for Ag and Cu, respectively. The calculations were considered to be converged with energy and force tolerances of 1 × 10⁻⁴ Ry and 1 × 10⁻³ Ry/a.u., respectively. The energy difference between decahedron (Dh) and icosahedron (Ih), defined as E_Dh − E_Ih, at the DFT/PBE level is +3.87 eV, +2.55 eV, and −2.56 eV for Cu, Ag, and Au, respectively. The corresponding values according to the Gupta potential are +1.57, +0.46, and −1.86 eV. Both DFT/PBE and Gupta, therefore, show the same trend: Ih is energetically preferred in Cu and Ag, while Dh is favored in Au. Based on these results, we believe that Gupta potentials are reliable for analyzing structural trends between Cu, Ag, and Au metal nanoclusters. The use of this model will allow a thorough sampling of the energy landscape, which would be hardly feasible with DFT. A detailed comparison of Gupta’s potential with DFT calculations is provided in the Results and Discussion section, which allows us to assess its performance and limitations.

Before the PTMD simulations, we calculate the global minimum at each size using a basin hopping Monte Carlo (BHMC)^19,37,79 optimization search. For each size, we run five independent search simulations with at least 2.5 × 10⁵ basin hopping steps.

The detailed procedure of the combined method of PTMD + HSA is described in a previous study.¹⁹ Here, we only recapitulate it briefly. In the PTMD simulations, there are two fundamental parameters: the number of replicas (M) and the temperature, T_m (m = 1, 2, 3, …, M) of each replica. All the replicas are in a canonical ensemble (NVT), and the exchange of configurations between a pair of replicas is attempted at specific intervals. The number of replicas is chosen such that we have at least 20–30% acceptance of the replica swaps. This is achieved by calculating an approximate caloric curve to identify the melting range and then adjusting the number of replicas and their temperatures to achieve the desired swap acceptance rate. All the PTMD simulations have been carried out in the Large-scale Atomic/Molecular Massively Parallel Simulator (LAMMPS).⁸⁰ We use a time step of 5 fs for the molecular dynamics evolution, and replica swaps are attempted every 250 ps. They are either accepted or rejected according to a Metropolis-like criterion.

We begin the PTMD simulations with all the replicas having the same structure, either a global minimum or a low energy structure. After discarding the initial phase of PTMD (∼0.5 µs), we sample configuration at 125 ps after a swap attempt for a total time of about 1–2 µs.

We will present the results of Cu and Ag nanoclusters. We note that the structural distribution of Au nanoclusters has been previously reported by us,¹⁹ and we use it here to make a comparison with Cu and Ag. In addition, we compare in detail the structures of Au, Cu, and Ag, which were not reported previously. To begin with, we discuss the finite-temperature structural distributions and then make a comparison to highlight the differences and similarities between Cu, Ag, and Au clusters. The melting points of all the metal nanoclusters in the current work are reported in Table I. We identify the melting point by first constructing the heat capacity (C_V) curve from PTMD simulations. The melting point is then calculated as the peak of the C_V curve.

Cu has a strong preference for icosahedral motifs as compared to Ag and Au.^2,3 The global minimums of Cu₉₀, Cu₁₄₇, and Cu₂₀₁ are shown in Fig. 1. The global minimums of Cu₉₀ and Cu₁₄₇ are both icosahedra, with Cu₉₀ having C_2v point group symmetry. However, with the EAM potential, the global minimum of Cu₉₀ was predicted to be an icosahedron with C_s symmetry.⁹ The best structure of Cu₂₀₁ is a decahedron with C_s point group symmetry.

In the case of Cu₉₀, the icosahedron (Ih) is the dominant motif at room temperature, with a very small amount of twins, decahedra (Dh), and mix structures [Fig. 1(a)]. The mix structures comprise several different geometric types. Predominantly, the mix structures consist of icosahedral-based geometries that either have amorphous regions or the entire structure adopts a configuration similar to the 92-atom chiral structure^9,83 with two missing atoms. The remaining mix structures consist of polydecahedra (p-Dh), which have more than one local fivefold axis.^19,37,84 With increasing temperature, the proportion of mix structures increases at the expense of Ih and peaks before melting at ∼ 600 K. Qualitatively, HSA predicts similar structural changes in Cu₉₀. The agreement between HSA and PTMD is good at room temperature, but thereafter there are quantitative discrepancies. The increase in mix structures is rather slow, according to HSA. For example, at 600 K, PTMD predicts 71.2% mix structures, while HSA predicts only 20.6%.

At size 147 [Fig. 1(b)], the icosahedron, which is the global minimum, dominates in the entire temperature range according to both PTMD and HSA. This indicates a high thermal stability of the icosahedral motif at this size. Moving on to Cu₂₀₁ [Fig. 1(c)], again, the global minimum structure, a decahedron, dominates at room temperature, and its proportion decreases steadily with temperature. Icosahedra compete with decahedra at higher temperatures, with the maximum proportion of Ih observed at 700 K just before melting. HSA, on the other hand, predicts a significantly higher amount of Ih at this temperature (77.2% vs 33.7%). Interestingly, fcc and twin structures are almost absent in Cu₉₀, Cu₁₄₇, and Cu₂₀₁ clusters in the entire temperature range.

The global minimum structures of Ag₉₀ and Ag₂₀₁ (Fig. 2) are decahedra, with both structures having C_s point group symmetry. The ideal icosahedron is the global minimum of Ag₁₄₇. These results are consistent with the previously reported global minima at these sizes for Ag clusters.^10,11,15

Ag₉₀ exhibits interesting structural changes [Fig. 2(a)]. From the HSA results, it is evident that the global minimum decahedron undergoes a partial transition to twin and mix structures with increasing temperature, leading to a combination of Dh + twin + mix structures at 250 K. Considering the PTMD results, the proportion of Dh, twins, and mix structures remains constant up to ∼ 450 K. This is a case of one-to-many solid–solid transition²⁷ where one geometrical motif, the Dh, is replaced by a coexistence of Dh, twins, and mix structures. Above 450 K, the proportion of mix structures increases at the expense of Dh and twins. The mix structures are a combination of polydecahedra³⁷ and distorted icosahedra with amorphous regions. The structural changes in Ag₁₄₇ and Ag₂₀₁ [Figs. 2(b) and 2(c)] are fairly straightforward. In both cases, the global minimum motif (Ih for 147 and Dh for 201) dominates in the entire temperature range, with other motifs nonexistent or in extremely small proportions.

We have recently¹⁹ reported the structural changes in Au nanoclusters, and hence, we will only summarize them briefly here (see the supplementary material, Fig. S2). The global minimum structures of Au₉₀, Au₁₄₇, and Au₂₀₁ are fcc, decahedron, and fcc (ideal truncated octahedron), respectively. At size 90, the global minimum motif, fcc, is dominant at lower temperatures and competes with twin and mix structures. With increasing temperature, fcc structures decrease along with an increase in mix structures. In the case of Au₁₄₇, the decahedron (global minimum) remains dominant up to melting, along with small amounts of twin and fcc structures. Above 400 K, Ih and mix structures begin to appear with mix structures dominating close to melting. In Au₂₀₁, there is a solid-solid transition from the fcc global minimum (a truncated octahedron) to a Dh at a low temperature around 200 K. Thereafter, the Dh dominates up to melting along with a small amount of twins (∼10%).

We stitch together HSA and PTMD results in order to get the structural changes in the entire temperature range in a single plot. Data for Au are taken from Ref. 19. Figure 3 compares the available results for Cu, Ag, and Au at all temperatures and sizes. HSA and PTMD are stitched together at 300 K. At the temperature where HSA and PTMD are joined, their structural distributions have excellent agreement except for small jumps in the case of Ag₉₀ and Au₉₀, where the trends are anyway consistent. This shows that our approach of combining HSA and PTMD is fairly robust and validated across various metal systems.

There are broadly three categories of structural changes that can be observed: type-(i) the global minimum remains the dominant motif up to melting, where amorphous takes over; type-(ii) solid–solid transitions occur, either completely or partially, well below melting temperature, resulting in an entirely different dominant motif; type-(iii) solid–solid transitions gradually occur leading to the co-existence of multiple motifs. The cases Cu₁₄₇, Ag₁₄₇, and Ag₂₀₁ fall into the first category, while Au₂₀₁ falls into the second category. All other cases fall into the third category, but with some differences. In Au₁₄₇ and Cu₂₀₁, the coexistence between motifs is present in a relatively narrow temperature range close to melting, whereas in all clusters of size 90, coexistence is already found at low temperatures.

The results show that ideal geometries corresponding to the “magic” sizes are not necessarily energetically preferred. Here we considered two “magic” sizes, 147 and 201. At size 201, the truncated octahedron has perfect geometry. However, only Au has a truncated octahedron as the global minimum, while the decahedron prevails for Cu and Ag. Even in Au, the global minimum transforms to Dh, which remains the dominant structure at finite temperatures. On the other hand, at size 147, which corresponds to a perfect icosahedron, both Cu and Ag have this structure as the global minimum. However, the decahedron is the global minimum of Au₁₄₇ with some icosahedra appearing only above 400 K. At size 90, all three systems have a different geometrical motif as the global minimum—Ih for Cu₉₀, Dh for Ag₉₀, and fcc for Au₉₀—which remains dominant (Cu, Au) or competes with other motifs (Ag). The structural distribution of Cu reinforces the strong preference for an icosahedral motif in Cu clusters.

We have, so far, discussed how the various geometrical motifs compete with temperature. In this section, we characterize the structural features of the various motifs.

Typical structures of Cu₉₀ are shown in Fig. 4, along with their energies relative to the global minimum. The icosahedron is the dominant motif of Cu₉₀, along with minor amounts of twin and Dh. The twin structures of Cu₉₀ typically have stacking faults [second structure in Fig. 4(a)]. At higher temperatures, icosahedra resembling the 92-atom incomplete Mackay icosahedron with C_3v point group symmetry are observed. These structures have two surface vacancies at various positions on the 92-atom cluster resulting in Cu₉₀ icosahedra. An example is shown in the third structure in Fig. 4(a). As the temperature increases further, some of these icosahedra undergo a twist and transform to mix structures resembling the 92-atom chiral geometry with tetrahedral T symmetry [the fourth structure in Fig. 4(a)]. The 92-atom chiral structure is the global minimum^9,83 of Cu₉₂ and has also been experimentally confirmed to have T symmetry from a comparison of the photoelectron spectra of Na and Cu clusters.^85–87 Again, the chiral-like Cu₉₀ clusters have two surface vacancies.

In the case of Ag₉₀ [Fig. 4(b)], along with the conventional decahedron (first structure), we find decahedra with either one (third structure) or two (second structure) hcp islands. When the two hcp islands are adjacent to each other, a local decahedral axis is formed at the intersection, which can be considered a polydecahedron (p-Dh)⁸⁴ that has more than one decahedral axis. The twin motif, which competes with Dh, consists of either a single hcp plane (the fourth structure) or stacking faults (the fifth structure). Moving on to Au₉₀ [Fig. 4(c)], the twins predominantly have a single hcp plane, differently from Cu₉₀ and Ag₉₀. In addition, Au₉₀ decahedra have deeper reentrant grooves [see the arrow in Fig. 4(c)] compared to decahedra of Cu and Ag. This is consistent with the general trend found in Ref. 2. The decahedron can undergo surface restructuring resulting in a mix structure [see the fourth structure in Fig. 4(c)]. Consider the four atoms (1, 2, 3, and 4) shown before (top) and after (bottom) restructuring. The atoms 2, 4 are pushed apart, and the atoms 1, 3 come closer leading to a {100} like arrangement.

At size 147, we observe a gradual change in the nature of icosahedra from Cu to Ag to Au. With increasing temperature, the perfect icosahedra become defective, initially with a single vertex vacancy [the second structure in Fig. 5(a)] and, at still higher temperatures, multiple vacancies (the third structure). These same vertex vacancies are also observed in Ag₁₄₇ icosahedra [the second and fifth structures in Fig. 5(b)]. However, along with the vertex vacancies, we also observe “rosette”^57,58 defects where the vertex atom protrudes to join the five nearest neighbors on the surface to form a six-atom ring. These are highlighted in blue for Ag₁₄₇ in Fig. 5(b), where either two or three “rosette” defects occur together. Icosahedra in Au₁₄₇, which appear mainly above 400 K, almost always have “rosette” defects, as shown in Fig. 5(c). Au₁₄₇ decahedra at higher temperatures exhibit deep reentrant grooves compared to the global minimum [second and third structures in Fig. (5c)]. The twins in Au₁₄₇ predominantly have single hcp planes, as shown in Fig. 5(c).

Finally, at size 201, all three systems have Dh as the dominant motif at room temperature. In Cu₂₀₁ and Ag₂₀₁, the various decahedra that are observed are all obtained by differing arrangements of nine additional atoms on a magic sized 192-atom Marks decahedron. The nine additional atoms are indicated in red [see Figs. 6(a) and 6(b)]. The twins in Cu₂₀₁ have a significant amount of hcp regions that are either completely hcp or consist of stacking faults. At higher temperatures, we observe icosahedra, which are incomplete 309-atom icosahedra. In Au₂₀₁, Dh is the dominant motif. In this case, the best Dh [second structure in Fig. 6(c)] is different from the typical decahedra observed in Cu and Ag, which are formed by adding nine atoms to the 192-atom decahedron. Instead, the best Dh of Au₂₀₁ is highly asymmetrical with deep reentrant grooves. However, at higher temperatures, we do observe Dh structures similar to those of Cu and Ag [the third structure in Fig. 6(c)].

In addition to the structures discussed earlier, we also observe structures that are not straightforward to categorize. We refer to these as mix structures, which occur in greater proportions at the smallest size of 90. The typical mix structures at size 90 are shown in Fig. 7. In a polydecahedron (p-Dh),⁸⁴ more than one decahedral axis is present within the same nanocluster. Examples of Cu₉₀ and Ag₉₀ p-Dh consisting of three decahedral axes are shown in the first image of Figs. 7(a) and 7(b). On the other hand, p-Dh is highly uncommon in Au₉₀. Another type of mix structure has an icosahedral region along with a disordered region. All three systems exhibit these structures [second image in Figs. 7(a) and 7(b), and first image in Fig. (7c)]. A third type of mix structure occurs when local icosahedral features are observed within fcc/twin [final image in Fig. 7(a)] or decahedron [final image in Figs. 7(b) and 7(c)]. This type of structure is mainly observed in Au and is less common in Cu and Ag clusters. The proportion of mix structures is significantly lower at the larger sizes of 147 and 201. We observe structures similar to those at size 90, with icosahedra mixed with disordered regions being more dominant. A detailed analysis of mix structures in Au clusters has been discussed previously.¹⁹ 

The structural distributions of Cu, Ag, and Au presented so far correspond to the Gupta potential, which does not account for the electronic interaction between atoms. In order to assess the performance of Gupta’s potential, we make a comparison with DFT calculations. We used PAW pseudopotentials with three types of exchange-correlation functionals—Perdew–Burke–Ernzerhof (PBE),⁷⁸ local-density approximation (LDA),⁸⁸ and PBE revised for solids (PBEsol).⁸⁹ 

We choose highly probable motifs (either two or more structures per metal per size) depending on the structural distribution. For instance, Ih and mix are the most dominant motifs in Cu₉₀. In the case of Ag₉₀, three motifs coexist—Dh, twin, and mix. Hence, we chose the lowest energy Ih, mix for Cu₉₀, and Dh, twin, mix for Ag₉₀. All the Cu and Ag structures used for DFT calculations are shown in Fig. 8. For a given combination of metal and size, we measure the energy difference of each structure with respect to the global minimum predicted by Gupta potential. These values are reported in Table II for Gupta potential, DFT/PBE, DFT/LDA, and DFT/PBEsol.

In the case of Cu₉₀ and Cu₁₄₇, Ih has lower energy, according to both Gupta and DFT. However, for Cu₉₀, Ih wins by only ∼ 0.08 eV in comparison to $>$1 eV for all three DFT calculations. On the other hand, Cu₁₄₇ has a very good quantitative agreement with DFT. For Cu₂₀₁, Gupta potential predicts Dh to have lower energy than Ih, in contrast to DFT. In the case of Ag₉₀, DFT favors twin in comparison to mix and Dh. According to DFT, the energetic ordering is E_twin $<$ E_Dh $<$ E_mix. Gupta potential, on the other hand, predicts Dh to have the lowest energy among the three. There is a good agreement between Gupta potential and DFT for Ag₁₄₇ and Ag₂₀₁.

Moving on to Au, the various structures used for DFT calculations are shown in Fig. 9. In the case of Au₉₀, we observe a lack of consistency among the various exchange-correlation functionals. There is a good agreement between Gupta potential, DFT/LDA, and DFT/PBEsol, with all three predicting a lower energy for fcc vs twin. However, DFT/PBE predicts twin to be the lowest energy structure. For Au₁₄₇, we considered all the motifs (other than amorphous) given their co-existence before the melting region. Au icosahedra typically have “Rosette” defects. Hence, we also considered the regular closed-shell 147-atom icosahedron and referred to it as Ih-reg in order to assess the competition between them. The energetic ordering according to Gupta potential is E_Ih-reg $>$ E_Ih $>$ E_mix $>$ E_fcc $>$ E_twin $>$ E_Dh. First, Ih-reg has higher energy than Ih according to Gupta potential and DFT calculations confirming that Au favors defective icosahedra consisting of “Rosette” defects. DFT/PBE predicts Ih to have lower energy than Dh, while Gupta potential, DFT/LDA, and DFT/PBEsol predict the opposite. When it comes to mix vs Dh, Gupta potential disagrees with DFT calculations, which predict mix to have lower energy than Dh. However, it is interesting to note that the mix structure is indeed a Dh with local rearrangement of a few atoms near one of the reentrant grooves (see bottom of Fig. 9). Hence, we believe that the Dh motif will also dominate at the DFT level, in agreement with Gupta’s results. Finally, for Au₂₀₁, both Gupta potential and DFT predict the same energetic ordering: E_twin $>$ E_Dh $>$ E_fcc. However, Gupta potential has a lower energy difference compared to DFT. As a result, we anticipate that the solid–solid transformation from fcc → Dh will be delayed and occur at a higher temperature than predicted by Gupta potential.

Overall, we observe the following trends. At size 147, Gupta potential performs fairly well, more so for Cu₁₄₇ and Ag₁₄₇, which exhibit excellent quantitative agreement between Gupta potential and DFT. In the case of Au₁₄₇, Gupta’s potential does a good job. First, it predicts that defective icosahedra are preferred with surface “rosettes.” Second, Ih has higher energy than Dh and mix, according to both Gupta and DFT. The only difference is mix, which is a distorted Dh with local rearrangement near the reentrant groove, is energetically preferred over Dh at the DFT level. At size 90, there is a qualitative agreement between Gupta potential and DFT for Cu but not for Ag and Au. In the case of Ag₉₀, twin is preferred at the DFT level, while Dh is preferred according to Gupta potential. In Au₉₀, there is internal disagreement among DFT exchange-correlation functionals. However, given the very small energy difference (absolute values are about 0.1 eV or lower), we expect a similar competition between fcc and twin as observed with Gupta potential. Finally, at size 201, both Ag and Au exhibit a qualitative agreement with DFT (although they underestimate the energy differences). In the case of Cu₂₀₁, Ih is preferred at the DFT level as opposed to Dh according to Gupta potential.

Finally, in order to understand how the embedded atom method (EAM) pair potential model performs in comparison to the Gupta potential, we calculated the energy differences using EAM potentials for Cu,⁹⁰ Ag,⁹⁰ and Au.⁹¹ The results are reported in the final column of Table II. In the cases of Cu, Gupta, and EAM, they exhibit similar performances. On the other hand, EAM seems to perform marginally better in the case of Ag. According to EAM, Ag₉₀ predicts twin to have lower energy compared to mix, and Dh is similar to DFT. In the case of Au, EAM performs poorly in comparison to Gupta. The major drawback with EAM is that it predicts Ih to be the lowest energy structure for Au₁₄₇.

Based on these results, we find that model potentials are still a good guide to selecting the main structural motifs and discussing general trends between metals, but in some cases they fail to select the lowest-energy motifs in agreement with DFT. We note, however, that there is a case, Au₉₀, where there is qualitative disagreement even between different types of DFT calculations. Moreover, in general, there are quantitative discrepancies between the different exchange-correlation functionals, which would make it difficult to assign precise temperature-dependent isomer probabilities even at the DFT level.

In this work, we have applied a computational framework that we proposed recently¹⁹ to study the size- and system-dependent structural distributions of Cu, Ag, and Au nanoclusters. In this method, we combine harmonic superposition approximation (HSA) and parallel tempering molecular dynamics (PTMD) in a complementary manner and calculate the structures of metal nanoclusters in the entire temperature range from 0 K to melting. We considered three cluster sizes, −90, 147, and 201, in the range 1–2 nm, of which 147 and 201 are “magic” sizes.

To begin with, “magic” sizes are not necessarily “magic” in that the global minimum is not always the ideal geometrical motif at that size. The perfect icosahedron and truncated octahedron are the ideal geometries at sizes 147 and 201, respectively. However, only in three out of six cases (Cu₁₄₇, Ag₁₄₇, and Au₂₀₁) does the global minimum correspond to the ideal geometrical structure. The global minima of Au₁₄₇, Cu₂₀₁, and Ag₂₀₁ are all Marks decahedra. At size 90, all three systems have a different global minimum: icosahedron for Cu₉₀, decahedron for Ag₉₀, and fcc for Au₉₀.

The structural changes in these systems can be categorized broadly into three groups: type-(i) global minimum is also the dominant motif at finite temperatures up to melting; type-(ii) solid–solid transformations lead to a completely different motif; type-(iii) solid–solid transitions lead to a co-existence of two or more motifs. The majority of the cases belong to the second and third groups, which include Cu₉₀, Cu₂₀₁, Ag₉₀, Au₉₀, Au₁₄₇, and Au₂₀₁.

Icosahedra are extremely dominant, with almost 100% abundance in Cu₁₄₇ and Ag₁₄₇ right up to melting. Similarly, decahedra are the dominant motif in Au₂₀₁ up to melting. In the cases of Cu₉₀ and Au₉₀, we find a significant proportion of mix structures close to melting. Although decahedra are dominant in Cu₂₀₁, we find significant amounts of icosahedra beyond 400 K. Finally, in Au₁₄₇, the proportion of Dh decreases gradually, and we find small amounts of twin, fcc, Ih, and mix structures co-existing at higher temperatures.

In contrast, Ag₉₀ and Au₂₀₁ undergo solid–solid transformations. Ag₉₀ exhibits a partial transformation Dh → Dh + twin + mix between 100 and 150 K. Beyond 150 K, the proportion of Dh, twin, and mix structures remains approximately constant up to 450 K, indicating the co-existence of multiple motifs. In the case of Au₂₀₁, fcc transforms to Dh below 200 K resulting in almost 100% Dh at room temperature, which remains dominant up to melting. In both instances, the solid–solid transformation occurs well below room temperature (⁠$<$200 K). As a result, it is non-trivial to predict the finite-temperature structures from the global minimum alone.

We also observed system specific differences across the three metals. Cu has a stronger preference for icosahedral structures. This is evident from the almost 100% abundance at sizes 90 and 147. While at size 201, a significant amount of icosahedra are observed above 400 K, which peaks around 700 K with ∼33%. In the case of Ag, icosahedra are mainly observed at the “magic” size of 147, where they occur with almost 100% abundance. Au, on the other hand, disfavors icosahedra, with icosahedra observed mainly at size 147 in small proportions beyond 400 K. Another interesting feature is the gradual change in nature of “rosette” defects in icosahedra at size 147 from Cu to Au. “Rosette” defects are completely absent in Cu but appear at higher temperatures in Ag. However, typically, almost all the icosahedra in Au have “rosette” defects. In contrast to Cu and Ag, decahedra in Au have deeper reentrant grooves.

Finally, a comparison of the performance of the Gupta potential with DFT reveals a few limitations of interatomic pair potentials. We observe a good agreement between Gupta potential and DFT at size 147. In other cases (Cu₉₀, Ag₉₀, Ag₂₀₁, and Au₂₀₁), the energetic ordering of the considered motifs is the same according to both Gupta potential and DFT, with Gupta energy differences being underestimated. In the case of Au₉₀, Gupta potential agrees with DFT/LDA and DFT/PBEsol but not with DFT/PBE. Finally, Gupta potential fares poorly in the case of Cu₂₀₁ since it predicts Dh to prevail over Ih. However, according to DFT, Ih should prevail over Dh. Notwithstanding these limitations, interatomic pair potentials remain indispensable since the wide exploration of the energy landscape of metal nanoclusters at the DFT level is simply not feasible. It is instructive to first obtain the structural distributions using interatomic pair potentials, e.g., Gupta, as performed in the current work, followed by DFT calculations to understand the limitations of the structural distributions. For instance, we observed that the Gupta potential predicts a complete solid–solid transformation from fcc → Dh below room temperature for Au₂₀₁. However, the energy difference between Dh and fcc is lower than predicted by DFT (0.0524 eV for Gupta potential vs $>$0.2 eV for DFT). Based on this information, it can be inferred that the transformation from fcc → Dh may occur at a higher temperature than predicted by the Gupta potential. A further check of another model potential, EAM, shows an overall performance of the same quality as the Gupta potential, with a better agreement with DFT for Ag₉₀ and a poorer performance for Au clusters.

Our method can be easily applied to any size or system for which reasonable models for atomic interactions are available. As a result, this method enables one to estimate the equilibrium proportion of various geometrical motifs as a function of temperature, which can then be used to compare with the experimentally obtained structural distribution.^67–69 This allows one to verify if the experimentally observed structures are in equilibrium or kinetically trapped metastable structures.

Supplementary material contains Parameters of HSA, PTMD; Structural distribution of Au nanoclusters.

This work has been supported by the project “Understanding and Tuning FRiction through nanOstructure Manipulation (UTFROM)” funded by MIUR Progetti di Ricerca di Rilevante Interesse Nazionale (PRIN) Bando 2017 (Grant No. 20178PZCB5). M.S. and A.G. acknowledge financial support from MIUR “Framework per l’Attrazione e il Rafforzamento delle Eccellenze per la Ricerca in Italia (FARE)” scheme, Grant SERENA No. R18XYKRW7J. R.F. acknowledges the Progetto di Eccellenza of the Physics Department of the University of Genoa for financial support and the International Research Network Nanoalloys of the CNRS for networking support. The authors acknowledge PRACE for awarding us access to Marconi100 at Cineca, Italy.

The authors have no conflicts to disclose.

Manoj Settem: Formal analysis (equal); Investigation (equal); Writing – original draft (equal); Writing – review & editing (equal). Cesare Roncaglia: Formal analysis (equal); Investigation (equal); Writing – original draft (equal); Writing – review & editing (equal). Riccardo Ferrando: Conceptualization (equal); Supervision (equal); Writing – original draft (equal); Writing – review & editing (equal). Alberto Giacomello: Conceptualization (equal); Funding acquisition (equal); Supervision (equal); Writing – original draft (equal); Writing – review & editing (equal).

The data that support the findings of this study are available from the corresponding authors upon reasonable request.