Structural and Energetic Analysis of Copper Clusters : MD Study of Cun ( n = 2-45 )

Simulações usando a dinâmica molecular foram efetuadas, considerando-se um potencial empírico para investigar geometrias, padrões de crescimentos, estabilidades de estruturas e energias para clusters de Cu n (n = 2-45). Os clusters estáveis otimizados foram calculados pelo rearranjo via processo de colisão. O presente procedimento apresenta-se como uma alternativa eficiente para a identificação do crescimento de clusters e como uma técnica de otimização. Foi verificado que os clusters de cobre preferem formar estruturas compactas tridimensionais em determinadas configurações enquanto os sistemas de tamanho médio apresentam simetria esférica. Além disso, também foram observadas correlações entre os arranjos atômicos e os números mágicos dos clusters. Particularmente, verificou-se que Cu 26 tem uma estabilidade equivalente ao sistema Cu 13 .


Introduction
Clusters are quite different from solid-state materials.They are aggregates of nanoscale size, with an intermediate state of matter between molecules and bulk.They also exhibit a range of unusual physical and chemical properties, such as structural, electronic, and thermodynamic.Metallic clusters have been the subjects of intense research.2][3][4][5][6] In this respect, the changes of cluster properties as a function of size, such as evolution from small to large clusters, is one of the most interesting issues. 6ystematic structural studies represent the starting point for understanding other general cluster properties.[7][8][9][10] Unfortunately, determination of equilibrium structures, and of atomic arrangements in TM clusters, still remains as a challenging task.Moreover, any experimental investigation and production of isolated microclusters are extremely difficult.
Computational studies provide helpful atomistic level simulations by using density functional/ab initio calculations [11][12][13][14] or any accurate empirical model potential energy functions (PEF) with efficient methods. 15For ab initio calculations, in spite of providing accurate results, it is difficult to determine the lowest-energy structure of large size cluster due to a prohibitive computational demanding.Usually, the first aim is to obtain the lowestenergy minimum of the PEF using global optimization techniques.Computational simulations for predicting cluster properties have been regarded as powerful tools relative to the experimental difficulties.For example, Genetic Algorithms (GA) 16,17 and basin hopping (BH) [18][19][20] have shown to be reasonably accurate and are widely used for inspecting the global minimum of various empirical PEFs. 21As alternative methods, minima hopping (MH) for complex molecular systems 22 and simulated annealing (SA) for closed-shell systems 23 have also been employed.Other methods used as complementary tools can be also proposed for identifying any global minima of hypersurfaces. 24oreover, theoretical approaches can supply a set of very simple formulas.There are several proposed empirical PEF's in literature for various systems 25 which can be used for predicting cluster properties.
Copper nanoclusters are interesting and important in the field of Physics and Chemistry of TMs and their alloys due to their useful applications in nanoscopic devices and catalysts.There exist various experimental works for free copper clusters [26][27][28][29] and density functional/ab initio calculations 12,[30][31][32][33][34][35] for copper microclusters.For larger Cu clusters several relevant studies with model PEF have been also reported [36][37][38][39][40][41] in literature.Doye and Wales 42 applied Sutton-Chen type potentials to determine the global minimum structures of metal clusters by using a Monte Carlo (MC) minimization approach.Using an empirical PEF Bayyari et al. 43 obtained stable structures of Ni, Cu, Pd and Pt microclusters via Molecular Dynamic (MD) simulations.The existing literature on Cu clusters has mainly focused on the structural and electronic properties.For example, in a recent work, Grigoryan et al. 39 used embedded-atom method to obtain a detailed description of copper clusters.Erkoç 44 has investigated the effect of radiation damage on copper clusters by performing MD simulation using empirical PEF.Such a potential model is applied in the present work for describing the interaction between copper atoms.Erkoç also studied Cu n (n = 3-55) clusters at room temperature (T = 300 K) with the same pair potential by using a MC technique. 40In reference 45 the structural properties of Cu 50 , Cu 100 and Cu 150 nanoparticles have been studied by using a modified version of diffusion MC method and by applying an empirical pair potential.The compact spherical shapes for stable structures of these nanoparticles have been reported.
In this paper, computational results are presented for isolated medium size copper clusters containing up to 45 atoms.This work follows a similar approach in regards to previous work on gold 46 clusters.MD simulations have been performed using an empirical PEF for copper. 40,47The goals of this work are to establish an efficient optimization method and to further understand the structural implications of this PEF by identifying the characteristic structural motifs associated with the stable minima of copper clusters.In particular, a possible geometrical packing phenomenon was studied for Cu 2 -Cu 45 sizes.In order to predict their structural and energetic properties, rearrangement collision processes 6,15,24,48,49 have been applied in the fusion regime.Similar growing up procedure has also been applied for silver clusters. 50The growing of the structures of copper clusters and also the magic numbers were characterized.The magic size indicates that any cluster with a certain size is more stable than its neighboring clusters against dissociation or fragmentation.It was found that the 13-and 26-atom clusters are particularly stable and also there are several other structures that are relatively stable.The rest of this work is organized as follows: next section presents theoretical background and other sections contain analyses of the results and some brief conclusions, respectively.

Computational Background
All calculations have been carried out using classical MD methods for investigating the structures of copper clusters through the Cu+Cu n-1 (n 45) collision.It is possible to compute the total interaction energy of a N-particle system from the sum of suitable effective-pair interactions. 25,41,47The effective-pair PEF used here is 40,47 (1) with the parameters 25 A 1 = 110.766008,A 2 = 46.1649783, 1 = 2.09045946, 2 = 1.49853083, 1 = 0.394142248, 2 = 0.20722507, D 21 = 0.436092895, and D 22 = 0.245082238 for copper.In these parameters, the energy is in eV and the distance is in Å. Erkoç has reviewed various potentials used in atomistic simulations. 25He reported the present empirical PEF for FCC metal microclusters of copper, silver and gold. 25In reference 47 it was pointed out that this PEF satisfies the bulk cohesive energy, and the bulk stability condition for Cu element.Moreover, the present PEF has been used to simulate copper nanowires. 51In the present work, based on equation (1) the MD was performed and the Runge-Kutta algorithm, of 5 th and 6 th order, is used as the numerical integration in all calculations.In the trajectory integration, Cartesian coordinates are used for time dependent positions and moments of the particles.All trajectories were checked during the integration to produce energy conservations of the order of 10 -10 in step size control of the microcanonical simulations.
The most stable isomers can be found by localizing the low-lying minima on the PEF.At the beginning of the atom-cluster collision the initial potential energy of the system is equal to the target Cu n-1 cluster energy.While this building-up procedure the colliding atom was sent from an asymptotic region.Formation of the new cluster is expected when the translation energy of the new projectile atom is released to the cluster in order to produce lower energy for the whole new cluster.The initial center of mass motion was kept constant during the interaction and the collision occurs around the center of mass of the system.To avoid fragmentation and scattering, all collisions are performed with low energies to keep particles together in the fusion regime.There is not any typical collision energy value because it can be change depending on sizes, sites, impact parameters, etc. Collision sites on the target cluster are also effective for these regimes.For example, when the colliding atom hits the target cluster on any open sites, it may easily construct a new structure.Moreover, the orientation of the target clusters is randomly represented by Euler angles. 52After the interaction starts, the collision energy of the projectile atom is distributed amongst the kinetic energies of all particles in the system through the rearrangement.The most stable orientation is determined by following each trajectory, set by checking the potential energy of the system at 100 steps up to the end of 2 million steps.The newly generated configuration that has the energy nearest the minima in each trajectory set is kept, and after 10 4 relaxation steps, it is minimized by removing kinetic energy step by step to determine the corresponding structure at 0 K.That is, the particles are kept under the force generated from the potential through this simple energy minimization procedure.Thermal quenching to reach 0 K means the minimization of the kinetic energy 5 i.e. the total energy is equal to just the potential energy values.After finding the new cluster, it is used in a new collision and this procedure is repeated in order to find new larger clusters.This procedure is repeated for randomly selected orientations of five initial configurations.As in our previous work 46,49 this methodology was applied to investigate the structures and the possible growing mechanism.

Results and Discussions
Optimized possible stable structures of copper clusters up to 45 atoms are presented in Figures 1 and 2. Most of the structures are similar to previously reported geometries of LJ clusters. 53,54The structures of small clusters up to Cu 26 are displayed in Figure 1.The well-known primitive geometries for 2-, 3-and 4-atom clusters are small enough to allow possible minima to be constructed directly.A regular tetrahedron is the most stable geometry of Cu 4 within T d symmetry.Its bond length and binding energy values have been calculated as 2.52 Å and 0.87 eV/atom, respectively.Such findings are found in similar studies using different empirical potentials.Grigoryan et al. 39 reported the same symmetry for Cu 4 via model potentials, an embedded atom-method (EAM) study.If one takes a DFT computations, such as Li et al., 55 one finds linear for Cu 3 and planar structures for Cu 4 and Cu 5 .In the present work, a trigonal, an octahedron, and a pentagonal bipyramids are predicted as ground state structures for Cu 5 (D 3h ), Cu 6 (O h ) and Cu 7 (D 5h ) clusters with 1.06, 1.26 and 1.38 eV/atom binding energies, respectively.In previous studies, 40 Erkoç and Shaltaf generated Cu clusters by using MC computation.
The structures, found here by using rearrangement collisions of clusters, are in agreement with their results 40 and also the results in reference 39.In this building-up procedure, from Cu 9 to Cu 12 clusters, the ground state geometries are in a growing pattern based on icosahedrons packing through filling the triangular open sites of 7-atom In a similar way, up to Cu 25 the clusters prefer to grow from the low coordination and more reactive sites on the equatorial region.A similar behavior was also observed in other studies of gold clusters as described in references 46 and 56.Cu 26 has an interesting view of crossed shape of 19-atom geometry.These structures are often based on the double icosahedrons geometry, with the additional atoms attached to various positions on 19-atom cluster as observed for titanium, vanadium and chromium clusters. 6Increasing the number of atoms on the surface of the cluster leads to some structural distortions of the basic building elements.This behavior was also verified in the previous studies of gold 46 and iron 49 clusters.As cluster size increases further, it becomes increasingly difficult to visualize the growth pattern.For microstructures consisting of a few atoms, it is easy to get a new structure by binding over a favorable open site.However, the new geometry for larger clusters may go in different local isomers of the new configuration.This is due to their dislocated structures and high symmetries. 46ost of the clusters, in general, might have various local minima corresponding to the absolute minimum energy of the PEF of a many-particle system.Metal clusters are extremely floppy, usually leading to numerous minima in the potential hypersurfaces.The considerable adsorption regions are atop, bridge, and hollow sites on cluster surfaces.From rearrangement structures, computationally it takes a long computational time to obtain the stable geometry of the clusters, unless the colliding atom hits the target at a suitable site.Therefore, the target position was randomly changed in order to improve the efficiency of the collision process.The cluster formation mechanism can be analyzed by investigating their preference for a growing pattern.As presented in Figure 2, from Cu 27 up to Cu 45 clusters, atoms prefer to fill favorable hollow sites (the most favorable adsorption site) on the target clusters.It is observed here that the new optimized structure grows from the hollow site of the previous smaller cluster.All configurations led to the migration of the colliding atom from the on-top or bridge sites to the hollow sites.Due to their reactivity, the adsorbing atoms were generally introduced onto low coordination copper atoms.As a consequence of this pattern, filling the hollow and more reactive sites, all low coordination points on the equatorial region of Cu 19 are covered one by one.Finally, a closed shell structure of Cu 34 is formed with 3 new five-atom rings, as can be observed in Figure 2. The structural evolution of the Cu 35 cluster in this formation pattern is demonstrated in the form of placing an atom on more reactive hollow site of the Cu 34 geometry.As a result, the new larger sizes will continue by filling the open sites of the surface of 34-atom cluster, such as the growing pattern of Cu 36 and larger clusters.With the increasing number of atoms the skeletal structure of Cu 34 gains stability and loses their original form inside of the larger clusters.The smaller sizes of the determined possible global minima have centered icosahedral morphologies.Octahedral, decahedral, and icosahedral morphologies have also been observed for the predicted low-lying structures, due to the increasing size.Larger sizes, up to Cu 45 , lead to more reactive, favorable hollow sites on the surface of clusters.Therefore, it becomes more complicated to determine the most stable structures.The DFT results [57][58][59] indicate that the sequence form sizes 34 to 45 atoms for copper clusters proceeds from the polyicosahedron of 34 atoms towards the anti-Mackay icosahedron of 45 atoms as determined the same structure in the present work.This common behaviour with DFT results indicates that the twobody potential used here can be considered as a reasonable PEF approach for this kind of analysis.
The calculated total energy values (E tot ) for Cu 2 -Cu 45 clusters are given in Table 1.The binding energies, the average interaction energy per atom in the cluster, versus the cluster size are plotted for the putative stable structures in Figure 3a.The decaying trend of the average binding energy with respect to the cluster size is an expected behavior for almost all metal clusters. 60s compared in Table 2, the calculated binding energies are, in magnitute, higher than the reported results in reference 40 and also closer to the experimental findings. 292][63] The average binding energy per atom in the cluster may be, therefore, expressed as a function of the cluster size, [64][65][66] (2) where the coefficients E v , E s , and E c correspond to the volume, surface, and curvature energies of the particles forming the cluster, respectively, and E e defines the energy origin. 67Fitting the data according to equation (2) gives 5.82, 11.45, 5.64 and 6.6 10 -5 eV/atom, respectively.In reference 40 this fit has been done with 3-terms for three parameters, E v = 1.736,E s = 2.727 and E c = 0.835 eV/atom.The experimental bulk cohesive energy is -3.49eV/atom 68 corresponding to calculated volume energy values.The differences between calculated and experimental values are 1.76 eV/atom in reference 40 and 2.33 eV/atom in the present calculation.However, the value determined here by our fitting procedure seems to converge to the experimental value faster than the previous work. 40That is, the corresponding size of the cluster for the experimental value here is about Cu 73 for equation (2) but it is in the limit of n in reference 40.Accordingly, the level-off value of binding energy here is 3.03 eV/atom for Cu 45 but in that work 40 it is 1.06 eV/atom for Cu 55 , far from the bulk cohesive value.A central issue in cluster physics is to identify particularly stable sizes.A detailed structural picture and the nonmonotonic variation in the cluster properties can be obtained by locating the global minimum as a function of size.Therefore, this can give information regarding the provided abundances of particularly stable clusters. 52Figures 3b and 3c show the first energy difference and the second finite difference (the stability function) of the total energy of the determined clusters (3) (4) as a function of the number of atoms, respectively.The peaks observed in Figure 3c correspond to the most stable structures (magic clusters) and the minima show the least stable sizes.Although it is known that the theoretical results of cluster stabilities are determined by the 2 E, this term has been assigned in the literature 14,15,19,39,46,49,67 as equivalent to the term of magic clusters as also indicated, for example, in references 6 and 40.In the actual work the same term is therefore used but one should remember that the correct is the second difference of the cohesive energy.The appearance of magic numbers for enhanced stability of the clusters and the fact that the clusters tend to form in spatial arrangements.The following magic numbers are observed: 13, 19, 23, 26, 29, 32, 34  In order to obtain further insight regarding the size dependence of structural growth and to make the magic sizes more deterministic, distributions of the atoms in their determined stable geometries have been investigated (Figure 4).Firstly, the radial distributions of the atoms are analyzed, which are displayed in Figure 4a.The radial distribution is the distance of each atom with respect to the center of mass of a Cu n cluster and it is given by, (5)   in which R i is the position of the i th atom.In the upper panel of Figure 4, all these distances are shown as a function of the cluster size.One aspect in the resulting diagram is the increasing radius of the clusters with increasing size.The largest distance to the origin (assumed as radius of the cluster) is increasing continuously with the increasing number of atoms.Some irregularities occur and in those cases the cluster radius decreases slightly by adding an atom e.g., the radius of Cu 5 is larger than Cu 6 , and Cu 10 has a larger radius than those of its larger neighbors Cu 11 and Cu 12 .The maxima in the largest distances correspond to the more reactive sizes.For example, the capped icosahedron form of Cu 14 , has more reactive sites due to the low coordination of the system.In particular, the trends of the radial distribution of the clusters (largest distances) have lower values identifying obviously for determined magic sizes.In most cases, this decrease is consistent with a reorganization of the system and an increase of the number of symmetry elements.In reference 69 Joswig and Springborg have noticed similar characteristics in aluminum clusters.Another aspect of the radial distribution analysis is that increasing the number of atoms per cluster leads to various different distances.It means that these clusters have lower symmetries than those with only a few different distances to the origin as stated in reference 69.A similar plot of the radial distribution for copper cluster was also observed in reference 39.It is also possible to identify more atomic shells using the radial distribution function.
The second shell of atoms is already established from Cu 13 , but for the smallest systems the inner shell contains just a single atom, which is placed very close to the center of the cluster.The microclusters up to the 13-atom cluster grow via the pushing of an atom to the center of the cluster.Cu 4 and Cu 6 have similar behavior since they are in regular tetrahedron and octahedron structures, respectively.In other words, all atoms are the same distance from the center.For magic sizes, the number of atoms of the inner shell can be easily observed due to their higher symmetric structures.For example, Cu properties, the mean displacements from the center of mass of the clusters are, as expected, slightly increasing due to the close packing phenomena.These results have also been observed in gold clusters. 46n Figure 4b, mean, minimum, and maximum values of the interatomic distances (pair displacement distributions) are demonstrated for Cu 2 -Cu 45 clusters as a function of the number of atoms.Maximum, minimum, and mean pair distances of atoms are the same, 2.52 Å, for the 4-atom cluster, due to its regular pyramidal geometry.The minimum pair distances decrease slightly while the mean pair distance increases with the increasing number of atoms.This is because the increase in the number of atoms leads to close packing of the system.When it reaches up to 45-atoms the mean and minimum values become 4.67 Å and 2.08 Å, respectively.On the other hand, the maximum pair distances have different trends in different size ranges.Structurally, different reorientations cause sudden increases and fluctuations in the maximum pair distances.For example, up to Cu 8 all structures are in different orientations.From Cu 8 to Cu 13 the icosahedrons packing based growing pattern based on pentagonal bipyramid structure of Cu 7 results in decrease for maximum pair distances.An addition of an atom to the triangular open sites of Cu 13 suddenly leads to a new increase in the maximum pair distance for Cu 14 , while the spherical structures of Cu 15 , with 6-atom rings, leads to a decrease in the maximum pair distance.Up to Cu 19 , the growing pattern is based on 13-atom geometry.The distance between two polar atoms of 19-atom cluster is the source of the rapid increase in the maximum pair distance from Cu 18 to Cu 19 .From Cu 19 up to Cu 26 one observes a slight decrease in the maximum pair distances.For the particular case of Cu 26 there is an interesting symmetric form, crossing shape of the two 19-atom clusters.After passing the Cu 26 structure there is a rapid increase in the maximum pair distance due to the new nonsymmetrical form of Cu 27 .There are slight fluctuations in regions 27-38 but no significant change in 37-45.In general, any changes are typically determined around the magic sizes.
Another alternative analysis used to investigate the growth mechanism in more detail was the calculation of density coefficients for number of atoms per volume based on the work in reference 42.For each cluster one can write (6)   that defines the proportionality relation for the cluster density, where n is the number of atoms and r n is the radius of the cluster corresponding to the largest value of the radial distributions in Figure 4a.Figures 5a and 5b illustrate the density coefficients and their stability functions i.e. the second finite difference of the calculated coefficients.Large fluctuations occur in the microcluster region.This means that they are reoriented through the changing of all atom positions.This fluctuation is smaller for medium size clusters because for larger clusters the orientation of the new clusters after rearrangement collision occurs on the surface atoms of the clusters.Inner structures of these clusters generally keep their previous geometries.The minima in the stability functions visually seem more symmetrical for close packing sizes, such as 4, 6, 13, 15, 23, 26, 29, and 34.There is an interesting result for Cu 15 that looks structurally like a magic size cluster, but it is not energetically favorable to form magic size cluster.This implies that spherical structures may not be energetically magic.In a similar way, Guo et al. 70 showed that a cluster with high symmetry is not always more stable than that with a lower symmetry.
Additionally, the moment of inertia (MoI) for these particular Cu 2 -Cu 45 clusters is analyzed in a similar way with regards to reference 71.The calculated results are presented in Figure 6.The values (I x , I y , I z ) of MoI with respect to the three components of the Cartesian coordinates are plotted as functions of the cluster size in Figure 6a.They have been calculated assuming the mass of the particles to be a normalized unit mass of 1. Therefore the units depend on the distance of the atoms in the clusters from the center of mass of the clusters.The equal values of three MoI show that cluster is in a spherical structure.As observed from the figure, the sizes n = 4, 6, 13 and 26 are spherical geometries.Second finite difference of the total MoI values are given in the stability graphs in Figure 6b.
The maxima in this figure demonstrate the relatively more spherical copper structures.The cluster sizes 13, 15, 23, 26, 29, 34 and 41 have magic size characteristics.Finally, mean values of the component dependent differences of MoI have been calculated by using absolute values of I x -I y , I y -I z and I z -I x .As presented in Figure 7, Cu 4 ,

Conclusions
In this paper a rearrangement collision procedure has been systematically performed to find likely global minima for the free Cu n clusters, in the size range of n = 2−45.These studies were based on the PEF proposed by the Erkoç. 40,47It has been shown, by using MD and energy minimization techniques, that copper clusters prefer to form three-dimensional compact structures and the five-fold symmetry appears on the spherical clusters.Particularly, all structural and energetic results show that the Cu 26 cluster has relatively more stable and spherical features.Therefore, it may be assumed to be another wellknown putative magic size.Additionally, high symmetry clusters are not always more stable than those with lower symmetries.As a result, the PEF can be used for qualitative structural analysis of medium size clusters such as for determinations of magic sizes.The rearrangement collision approach may be used as an alternative procedure for the investigation of possible structures of atomic clusters.The procedure can easily be applied to other cluster systems with different interatomic PEFs.In addition, the selected calculations of moment of inertia, radial, pair, and density coefficient distributions can be considered as efficient tools for structural analysis.

Figure 3 .
Figure 3. a) binding energies b) first and c) second finite difference of the total energies for the predicted stable copper clusters.

Figure 4 .
Figure 4. a) radial distributions and b) interatomic distances of atoms for the predicted stable copper clusters.

Figure 5 .
Figure 5. a) density coefficients and b) their stability values.

Figure 6 .
Figure 6.a) the moment of inertia (MoI) with respect to x, y and z components and b) the stability function of the total MoI values for copper clusters.

Cu 6 ,
Cu 13 and Cu 26 are exactly in spherically symmetric geometries.

Table 1 .
The calculated total energy values (E tot ) for Cu n (n = 2-45) clusters

Table 2 .
Comparision of the calculated binding energies (eV/atom) with previous available theoretical and experimental values for Cu n (n = 2-8) clusters 19and Cu 26 have mainly four distances.There is a dominant deviation of the central atom from Cu 12 to Cu 26 clusters.It reaches the biggest value in this medium size region.Cu 26 is the turning point for the central atom because the growing structure at this point has half filled equatorial sites of the Cu 19 cluster.Even though the closest and the largest distances have absolutely different