Estimating the Impact of an All-Electron Basis Set and Scalar Relativistic Effects on the Structure , Stability , and Reactivity of Small Copper Clusters

Basis sets of valence double and quadruple zeta qualities and the Douglas-Kroll-Hess (DKH) approximation are used to estimate the impact of an all-electron basis set and scalar relativistic effects on the structure, stability, and electronic properties of small neutral copper clusters (Cun, n ≤ 8). At the Becke three-parameter for exchange and Perdew-Wang 91 for correlation (B3PW91) non-relativistic and relativistic levels of theory, the bond length, binding energy, ionization potential, electron affinity, chemical potential, chemical hardness, and electrophilicity index are calculated. The results show that the agreement with experiment improves significantly when the DKH Hamiltonian combined with an all-electron relativistic basis set is used. Polarizabilities and hyperpolarizability are also reported. At the B3PW91 level, all-electron basis sets are shown to be more reliable than effective core potential valence basis sets in the determination of the second hyperpolarizability of copper clusters.

The goal of most studies has been to examine how the properties of a cluster evolve with size.These properties include the geometric structures, binding energies, ionization energies, and the highest occupied molecular orbital (HOMO) and lowest unoccupied molecular orbital (LUMO) energies.
It is not a trivial task to calculate properties of large clusters because the computational cost increases rapidly with increasing cluster size.Several computational strategies have been tested.Among them, we call attention to the effective core potential (ECP) approach. 15An appropriate ECP with suitable valence basis sets reduces the computational effort because only valence electrons are explicitly treated.
We recall that the ECP approach in combination with the LANL2DZ valence basis set have, to date, been used in most calculations of the structures and stabilities of small copper clusters. 10,13However, relativistic effects have been omitted in these calculations.
To obtain reliable theoretical results for copper clusters, both the relativistic and electronic correlation contributions must be considered for high-quality all-electron basis sets.
The scalar relativistic correction can be substantial, even for molecules that contain first-row elements. 16he Douglas-Kroll-Hess (DKH) approach [17][18][19] has been successfully used to estimate such scalar relativistic effects.
In this paper, we concentrate on the treatment of molecular scalar relativistic effects using the DKH method.We employ the Becke three-parameter for exchange and Perdew-Wang 91 for correlation (B3PW91) 20,21 hybrid functional with all-electron segmented contraction non-relativistic basis sets of valence double and quadruple zeta qualities plus polarization functions (DZP and QZP), 22,23 which are also available in the DKH recontraction, 24 and the augmented DZP and QZP (ADZP and AQZP) 22,25 sets for the Cu element.At the non-relativistic and relativistic levels of theory, the bond length, binding energy, HOMO and LUMO energies, ionization potential, electron affinity, chemical potential, chemical hardness, electrophilicity index, and (hyper)polarizabilities of the ground state of neutral copper clusters are calculated and compared with the experimental and/or theoretical data reported in the literature.

Methodology
All calculations have been performed with the Gaussian 09 code. 26hroughout the calculations, the DZP, ADZP, QZP, and AQZP and DZP-DKH, QZP-DKH, and AQZP-DKH basis sets are used with the non-relativistic and relativistic second-order DKH (DKH2) 27 Hamiltonians, respectively.These sets are available at the web site, 28 where basis sets are provided in properly formatted forms for the commonly used molecular program packages.
Initially, at the B3PW91/DZP and DKH2-B3PW91/ DZP-DKH levels of theory, the equilibrium geometries of a given cluster are determined.From these optimized geometries, the other properties are then calculated using the QZ basis sets.

Results and Discussion
From equations 3 and 4, it is clear that the IP and EA are fundamental to obtaining a reliable estimation of chemical potential and hardness.These properties are sensitive to the functional used in the DFT calculations.In light of this information and because the B3PW91 functional has been successfully applied in nonrelativistic calculations of the IP and EA of atomic and molecular systems, 10,13 this functional seems to be a natural choice.

Structures of copper clusters
Table 1 displays the bond lengths obtained from different levels of theory for ground state copper clusters of up to 8 atoms.The spin multiplicities are singlets for the even-numbered clusters and doublets for the odd-numbered clusters.The optimized ground state structures show that Cu 3 -Cu 6 and Cu 7 -Cu 8 are planar (2D) and 3D (see Figure 1), respectively.
]10,13 The B3PW91/LANL2DZ bond length is 2.254 Å, 3 whereas the B3PW91/DZP computation has a bond distance of 2.250 Å.The bond distances obtained with the two non-relativistic calculations overestimate the experimental value of 2.2197 Å 34 by 0.03 Å.This discrepancy is due to the neglect of scalar relativistic effects: for the close-shell copper dimer, the first-order spin-orbit effect is zero because the total electronic orbital (L) and spin (S) angular momenta are zero, and consequently, the total electronic angular momentum, J = L + S = 0.In contrast, our relativistic calculation for the Cu 2 ground state shows a contraction of the bond distance by 0.031 Å.Thus, the inclusion of relativistic effects brings the DKH2-B3PW91/DZP-DKH bond distance into excellent agreement with the experimental one.For Cu 2 , it should also be noted that the B3PW91/LANL2DZ harmonic frequency 3 is equal to 260 cm -1 and that the accordance between the theoretical and experimental (264.55 cm -1 ) 34 harmonic vibrational frequencies improve from B3PW91/QZP (254.82 cm -1 ) to DKH2-B3PW91/QZP-DKH (269.87 cm -1 ).From these results and our experience with relativistic calculations, a pragmatic strategy for the calculations of copper clusters consist of the geometry Vol. 27, No. 1, 2016   optimization at the DZ level followed by single point energy and (hyper)polarizabilities calculations at the QZ level.
From Table 1, it is clear that the bond length always decreases going from the LANL2DZ valence electron basis set to the DZP all-electron basis set.This decrease is small for the dimer (0.004 Å), but it increases for the other clusters and achieves 0.16 Å for the trimer.However, a systematic reduction by approximately 0.03 Å is observed when relativistic effects are included.These results show the importance of using scalar relativistic effects along with an all-electron basis set to obtain reliable copper cluster structures.
From Cu 2 to Cu 8 , the DKH2-B3PW91/DZP-DKH average Cu-Cu bond lengths (<r Cu-Cu >) are 2.219, 2.347, 2.343, 2.352, 2.354, 2.414, and 2.440 Å, respectively.The variation of <r Cu-Cu > with n shows clearly structural transitions occurring from one-to two-dimension and from planar to three-dimension clusters.In addition, it is quite evident that the average bond length for the cluster size increases consistently and approaches the experimental distance in the bulk metal (2.556 Å). 35 To achieve that distance, it will be necessary to consider larger clusters.

Binding energy
From the optimized structures of the copper clusters, the BE values are calculated at the non-relativistic and relativistic levels.The values are listed in Table 2.The B3PW91/LANL2DZ 10 results are also included.At any level of theory, the BE increases with the size of the cluster.As one can see from Table 2, the BE decreases from B3PW91/LANL2DZ to B3PW91/QZP.The difference between the corresponding non-relativistic results increases    36 and Knudsen cell mass spectrometry (3.04 eV) 37 are available.Table 2 shows that our DKH2-B3PW91/QZP-DKH binding energies are each in satisfactory and excellent agreements with these direct measurements.
For larger clusters, the BE values were obtained from collision-induced dissociation (CID) experiments of anionic 7 and cationic 8 copper clusters.These energies are also included in Table 2.The BE values derived from the CID of the cationic clusters are considerably smaller than those of the anionic clusters.Table 2 shows that our non-relativistic and relativistic BE values agree well with those derived from the CID of cationic and anionic copper clusters, respectively, and that they are within the uncertainty bars of the experiment.As expected, the agreement between theory and experiment improves going from B3PW91/LANL2DZ to DKH2-B3PW91/QZP-DKH.
From Cu 2 to Cu 8 , the DKH2-B3PW91/QZP-DKH binding energies per atom (BE/n) increases monotonically with the cluster size: 0.951, 1.006, 1.327, 1.447, 1.629, 1.733, and 1.817 eV.In the limit of n going to infinity, it is expected that BE/n reaches the metal cohesive energy (3.50 eV for copper). 35However, the clusters studied in this work are still too small to yield a satisfactory estimate of the copper cohesive energy.
The HOMO-LUMO gap of the transition metal cluster is an important qualitative characteristic for studying the modification in band structure according to the cluster size.The HOMO-LUMO gap is also used to determine the capacity of a molecule or cluster to participate in chemical reactions.The value of DKH2-B3PW91/ QZP-DKH HOMO-LUMO gap has been calculated for all cluster sizes, and the values are 2.510 (Cu), 3.459 (Cu 2 ), 1.269 (Cu 3 ), 2.094 (Cu 4 ), 1.632 (Cu 5 ), 3.314 (Cu 6 ), 1.479 (Cu 7 ), 2.622 (Cu 8 ) eV.One can easily note that the value of HOMO-LUMO gap changes abruptly with cluster size and that clusters with an even number of atoms have a larger HOMO-LUMO gap, which indicates that they must be less reactive than clusters with an odd n.The extra stability exhibited by even copper clusters is due to their closed-shell configuration.

Ionization potential and electron affinity
The experimental [2][3][4][5] vertical IP and EA and the corresponding theoretical values obtained from equations 1 and 2 are presented in Table 4.The former are given with their error bounds.The general trend of the IP shows a notable oscillation.Clusters with an odd number of atoms have an IP smaller than do those with an even n.This leads to maximum and minimum values for even-and odd-numbered clusters, respectively.The EAs also show a characteristic oscillation, but unlike the IP, the minimum and maximum occur for the even-and oddnumbered clusters.
The explanation for these oscillatory trends is that even-and odd-numbered copper clusters have closed and open shells, respectively.The latter systems have only one electron in the HOMO, which is easier to remove than an electron from the HOMO (doubly occupied) of closed-shell systems.Consequently, the ionization potentials in odd-and even-numbered systems display minimum and maximum values, respectively.Exactly the opposite occurs for electron affinities.Because the odd numbered copper clusters can more easily acquire an electron in the openshell HOMO than in the LUMO of closed-shell systems, the maximum and minimum electron affinities occur for odd-and even-numbered systems, respectively.Similar to the HOMO and LUMO energies, the IP increases from B3PW91/LANL2DZ to B3PW91/QZP and from B3PW91/QZP to DKH2-B3PW91/QZP-DKH.A similar trend is observed for the EA.The DKH2-B3PW91/ QZP-DKH ionization potentials gave the smallest errors and, in general, are within the experimental uncertainty bars.Compared with the EA experimental data, the results obtained from all-electron basis sets significantly improve the agreement.However, the DKH2-B3PW91/AQZP-DKH electron affinities are more accurate.To the best of our knowledge, these are the best DFT results reported to date in the literature. 10,13Even so, for Cu 3 and Cu 6 , the accordance between theory and experiment is still poor.It should be noted that the difference between AQZP-DKH and LANL2DZ electron affinities always exceed 22% (for Cu 6 ) and can be as large as 79% (for Cu).These results confirm the importance of accounting for relativistic scalar effects with an all-electron relativistic basis set in calculations for metal clusters.

Chemical potential, molecular hardness and electrophilicity
In Table 5, experimental and theoretical data of m, η, and ω obtained from equations 3, 4, and 5 are displayed.The error bars in the experimental data are due to the errors resulting from the use of the experimental IP and EA.
The chemical potential is related to charge transfer from a system to another with a lower value of m.Thus, it is expected that the odd clusters present maximum values of m because they have an open shell and that after the transfer of one electron, they will close their electronic shell and will be thus more stable than the original open-shell clusters.In Table 5, the theoretical and experimental data display this oscillating behavior, with a local maximum for clusters with an odd n.Cu 3 has a high value for m because it is able to transfer an electron to close its electronic shell.The B3PW91/LANL2DZ chemical potentials are underestimated, and the DKH2-B3PW91/QZP-DKH/ AQZP-DKH procedure gives the best agreement with the experiment data.Except for Cu 3 , the latter results are within the uncertainty bars of the experimental data.In many cases, chemical hardness may be used to characterize the relative stability of molecules and aggregates.Because the principle of maximum hardness (PMH) 38 declares that molecular systems at equilibrium present the highest value of hardness, it is expected that the hardness displays an oscillating behavior with local maxima at the even clusters.The theoretical and experimental hardness values shown in Table 5 indicate that the required even-odd oscillating features with stable clusters (n even) are harder than their next door neighbors (n odd).In general, the three theoretical approaches give similar hardness results.
Electrophilicity measures the energy stabilization when the cluster acquires an additional electronic charge from the surroundings.From Tables 4 and 5, it is evident that ω is closely related to EA and exhibits similar variation in changing cluster size.Stable systems, i.e., less reactive systems, are expected to have low electrophilicity values because they are less likely to acquire additional electronic charge to become more stable.Along the theoretical and experimental series, an even-odd oscillatory trend is observed, and the odd copper clusters present a maximum electrophilic value because their electronic shells are closed when they receive an electron.The only exception occurs for Cu 3 because the inaccuracy of our computed EA values has a large effect on the corresponding electrophilicity results.Odd copper clusters are more inclined to accept an electronic charge from the environment than are even clusters.
Except for Cu 2 at the DKH2-B3PW91/QZP-DKH/ AQZP-DKH level, the theoretical ω values are underestimated, but the agreement with the experimental data always improves when the all-electron basis set and the scalar relativistic correction are considered.With exception of Cu 3 , the relativistic results are again within the experimental uncertainty bars.We believe that it is the first time that such accuracy has been achieved compared with previous theoretical and experimental data. 10The difference between corresponding B3PW91/LANL2DZ a n d D K H 2 -B 3 P W 9 1 / Q Z P -D K H / AQ Z P -D K H electrophilicity values is large (ca.8.3%).

(Hyper)polarizability
Static polarizability of clusters is an important property because it is proportional to the number of electrons of the system and because it is sensitive to the structure and shape of the system.
The mean dipole polarizabilities, polarizability anisotropies, and second hyperpolarizabilities of copper clusters up to the octamer reported in this work are calculated at the B3PW91/DZP geometries (see Table 1).Thea, Da, andg results are collected in Table 6.At any level of theory, one can verify that when going from the atom to the octamer, the mean dipole polarizabilities of copper clusters increase monotonically and present the expected proportionality with n.In contrast, a/n oscillates when going from Cu to Cu 4 and decreases from Cu 5 .The polarizability anisotropy for Cu n (n ≤ 8) increases as we move from the dimer to the hexamer and decreases as we move to the heptamer and octamer, i.e., a maximum value is found at the hexamer.Exception occurs at the B3PW91/LANL2DZ level, where the pentamer gives a maximum value.These data, in conjunction with the topologies of the copper clusters (see Figure 1), show that the polarizability anisotropy is directly related to the cluster structure.For example, at the planar clusters, it increases with the number of copper atoms.It is important to note that as the cluster structures become compact, as occurs in the cases of the heptamer and octamer, the polarizability anisotropy values decrease.The value for the octamer is slightly larger than the value of the dimer, which has an open structure.This tendency is similar to that reported by de Souza and Jorge, 39 who studied static polarizabilities on lithium and sodium clusters and verified that the polarizability anisotropies achieve minimum values for clusters containing 2 and 8 atoms and maximum values for planar ones.Except for Cu, the LANL2DZ and AQZP mean dipole polarizabilities are very similar, whereas the opposite occurs when we compare the corresponding BP86/TZVP-FIP1 (all-electron basis set of triple zeta valence quality augmented with seven field-induced polarization functions) 6 and B3PW91/AQZP results.In this case, the difference is 11 atomic units (a.u.) for Cu 3 , but this discrepancy reduces to less than 2.9% with the cluster size enlargement, showing a smaller dependence with the procedure used.However, with a few exceptions, the Da values are very sensitive to the basis set.This finding demonstrates that there is a strong dependence between basis set and anisotropy.In this case, a deeper analysis is necessary.
Neogrady et al. 40 computed the polarizability of a few metals using high-level-correlated-relativistic calculations.For the copper atom, they found 46.50 a.u., which is in excellent agreement with the B3PW91/AQZP result.The LANL2DZ result 13 is overestimate by 7.3%, but more recently, Calaminici et al. 6 reported a value of a = 47.02 a.u., which is 1.1% higher than the value reported by Neogrady et al. 40 The MP2/AVDZ4 all-electrona (161.24a.u.) and Da (138.68 a.u.) results for Cu 4 reported by Maroulis and Haskopoulos 41 agree quite well with B3PW91/ AQZP.This reinforces the idea that the all-electron AQZP basis set is a reliable choice to carry out polarizability calculations of metal clusters.It is interesting to note that the CAM-B3LYP/AQZPa and Da values for Cu n (n ≤ 4) clusters that were reported recently by Martins et al. 25 are very similar to those computed in this work.
In Table 6, we give the B3PW91/AQZP second hyperpolarizabilities for the copper clusters from Cu to Cu 8 .One can observe that from Cu 4 , there is a significant increase in the hyperpolarizability values.For copper clusters, the results for this property in the literature are scarce.Based in observations reported by Shigemoto et al., 42 Maroulis 43 discussed the possibility of a negative second hyperpolarizability for Cu 2 .Using a high level of theory [CCSD(T)] along with the all-electron basis set [7s6p6d2f]), the author obtained 10 -3 ×g = 86 a.u., 43 which is in excellent accordance with the B3PW91/AQZP result of 87.30 a.u.For Cu 4 , the effect of the all-electron basis set on the hyperpolarizability was estimated previously. 44At the MP2 level, the Stevens/Basch/ Krauss ECP in conjunction with the modified triple-split basis set CEP-121G result (123.5 a.u.) 44 showed a significant discrepancy with the all-electron [5s4p4d] basis set value of 271 a.u. 41Compared with the allelectron CCSD(T)/[5s4p4d] second hyperpolarizability value of 10 -3 ×g = 293 a.u., 41 the B3PW91/AQZP value (282.41 a.u.) agrees very well, whereas the B3PW91/ modified CEP-121G is underestimated (230.4 a.u.). 44For the other clusters displayed in Table 6, we are not aware of any previous second hyperpolarizability results reported by other authors.

Conclusions
In this paper, we have calculated and rationalized the all-electron basis set and scalar relativistic effects on the structure, stability, and electronic properties of small copper clusters (Cu n , n ≤ 8).
At the B3PW91 level of theory, we verify that the bond length always decreases going from LANL2DZ (valence electron basis set) to DZP (all-electron basis set) and from DZP to DZP-DKH (relativistic all-electron basis set).The opposite occurs for the HOMO and LUMO energies, vertical ionization potentials and electron affinities, chemical potential, and electrophilicity.The binding energy and chemical hardness decrease as we move from LANL2DZ to QZP and increase as we go to QZP-DKH.
For the studied clusters, the all-electron basis set and scalar relativistic effects on the bond distance, LUMO energy, electron affinity, and electrophilicity are large.In these cases, we believe that our results are the most accurate reported to date in the literature.
From the B3PW91/QZP-DKH results, the following specific conclusions can be drawn: (i) except for a few cases, the theoretical values are within the uncertainty bars of the experiment; (ii) the binding energy per atom of the cluster increases with cluster size; (iii) the ionization potentials and electron affinities show oscillatory behavior for even-and odd-numbered clusters because of their closed-and open-shell HOMO; (iv) the electron affinity and electrophilicity show a similar variation with cluster size; (v) chemical potential and hardness present an even-odd oscillating feature in agreement with the chemical insight.
Finally, we verified that the polarizability anisotropy is more sensitive to the basis set than the mean dipole polarizability and that for Cu 2 and Cu 4 , the hyperpolarizabilities reported in this work are in good accordance with results obtained at the CCSD(T) level of theory.Thus, the B3PW91/AQZP procedure seems to be a good option on (hyper)polarizability calculations of metal clusters.

Table 1 .
Optimized bond lengths r ij (see Figure1) for Cu n clusters.The symmetry is given in parentheses 34From reference 10; b this work, all-electron basis set from reference 22; c this work, all-electron basis set from reference 24; d experimental value 2.2197 Å.34

Table 2 .
Binding energies of the fully optimized structures of copper clusters Cu n Estimating the Impact of an All-Electron Basis Set and Scalar Relativistic Effects J. Braz.Chem.Soc.130 from 0.145 eV for Cu 2 to 0.702 eV for Cu 8 .In contrast, the scalar relativistic effect increases the magnitude of the BE from 0.123 eV for Cu 2 to 0.913 eV for Cu 8 .For Cu 2 and Cu 3 , direct measurements of the atomization energies by fluorescence (2.08 eV) aB3PW91 results from reference 10; b this work (B3PW91), all-electron basis set from reference 22; c this work (DKH2-B3PW91), all-electron basis set from reference 24; d reference 7; e reference 8.

Table 3 .
HOMO and LUMO energies of the fully optimized structures of copper clusters Cu n aB3PW91 results from reference 10; b this work (B3PW91), all-electron basis set from reference 22; c this work (DKH2-B3PW91), all-electron basis set from reference 24.Vol.27,No. 1, 2016

Table 4 .
Vertical ionization potentials (IP) and electron affinities (EA) of the fully optimized structures of copper clusters Cu n