Dispersion Interactions in Small Zinc, Cadmium, and Mercury Clusters

Article pubs.acs.org/JPCA

Dispersion Interactions in Small Zinc, Cadmium, and Mercury Clusters Richard Hatz, Vesa Han̈ ninen, and Lauri Halonen* Laboratory of Physical Chemistry, Department of Chemistry, University of Helsinki, P.O. Box 55 (A.I. Virtasen aukio 1), FIN-00014 Helsinki, Finland S Supporting Information *

ABSTRACT: Dispersion interactions in small clusters of group XII (Zn, Cd, Hg) metal atoms are studied at the CCSD(T) level with triple-ζ basis sets. A pair potential model together with a least-squares fit to the interaction potential energy surface is used to calculate interatomic dispersion coefficients, which are found to be in good agreement with atomistic calculations. The angular dependence of the dispersion interaction, extracted by explicitly accounting for the leading order electrostatic and induction terms, is determined with the aid of a coupled spherical harmonic expansion. The distance dependence of the orientation averaged dispersion energy is modeled by an integral average of our potential model, which is able to robustly account for the average dispersion energy to a remarkable degree.

■

INTRODUCTION The long-range van der Waals (vdW) dispersion interactions are known to play a decisive role in a wide variety of fields.1−4 Although important, the vdW forces can also be difficult to model accurately. The most basic electronic structure calculation methods, the Hartree−Fock theory and some older functionals in density functional theory (DFT), are even known to completely lack the description of dispersion forces.5−7 Because of this, the dispersion forces have received a great deal of attention, and therefore several methods have been developed to accurately account for them. In DFT, for example, many different corrections, new functionals, and other schemes have been proposed to correct this shortcoming.8−10 Also among the pure ab initio methods, many different theoretical frameworks such as the symmetry-adapted perturbation theory11,12 (SAPT) and the linear response theory13,14 have been developed and successfully applied to the calculation of dispersion forces. However, even with the most modern tools available in computational chemistry, the accurate modeling of dispersion interactions can present difficulties, especially in cases where complicated or extended systems such as protein folding or molecular adsorption are under study.15−17 Besides a perturbative approach, the dispersion energy contributions to a molecular system can also be calculated by using a supermolecular method, where the interaction energy is calculated by subtracting the monomer energies from the total energy. The interaction energy is then further divided into distinct contributions such as long-range electrostatic, induction, and dispersion interactions and short-range interactions arising from the overlap of the wave functions between the interacting species. Because of its simplicity, this method has been extensively applied in many different systems,18−22 but the drawback is that the dispersion energy coefficients obtained are obviously dependent on the particular system, thus limiting © 2014 American Chemical Society

the applicability of this scheme in a more general context. The importance of being able to describe the dispersion energy contribution in a system-independent manner has been widely recognized, and many different kinds of methods have been applied to deal with this problem ranging from simple empirical formulas19,23 to rather sophisticated computational methods, e.g., in the DFT-D3 model by Grimme24 and co-workers and in the method proposed by Tkatchenko and Scheffler.25−27 In this study, which follows along the lines of two previous publications,28,29 we investigate one such simple scheme, where atomic dispersion coefficients are extracted from molecular interactions, and which could potentially be used to describe dispersion energy in many different kinds of systems. The promising results obtained previously for coinage metal clusters are now extended to group XII (Zn, Cd, Hg) clusters, with some interesting notions concerning the convergence of the dispersion energy fit to the long-range interaction energy potential energy surface. The results obtained here could, for example, be used to model the dispersion forces in larger group XII clusters, which, in light of previous work done in the field,30−35 constitute an interesting test case to study, e.g., the change in bonding from van der Waals-like to metallic.

■

THEORY Long-Range Interactions. Considering the behavior of the intermolecular interaction energy with respect to the distance R between the interacting species, the contributions can be divided into two distinct parts.36,37 The short-range contributions, arising from the wave function overlap, decay exponentially with respect to the distance and can be ignored at long distances. Received: May 22, 2014 Revised: July 7, 2014 Published: July 8, 2014 5734

dx.doi.org/10.1021/jp505023g | J. Phys. Chem. A 2014, 118, 5734−5740

The Journal of Physical Chemistry A

Article

For example, a molecule in the D∞h point group can only have one independent component of a multipole moment of any even rank.40 In such a molecule, the few lowest order contributions to the electrostatic energy are the quadrupole− quadrupole (EΘΘ), quadrupole−hexadecapole (EΘΦ), and hexadecapole−hexadecapole (EΦΦ) terms, described by the quadrupole moment Θ and hexadecapole moment Φ. The lowest order contribution to the induction energy is the quadrupole induced dipole interaction EΘα described by Θ and the dipole−dipole polarizability tensor α, which, for a molecule in the D∞h point group can have only two nonzero components, which are either parallel or perpendicular to the C∞ rotation axis, denoted α∥ and α⊥, respectively.38

This means that the interaction energy of the whole system at long intermolecular distances can be calculated from the properties of the interacting species alone. This energy can be further divided according to physical origin into a sum of the electrostatic, induction, and dispersion energies. With the help of the multipole expansion, these terms can be represented as a power series in R. Vlong =

∑ CiR−i

(1)

i

where the parameters Ci are functions of the geometry of the system and the multipole moments and polarizabilities of the interacting molecules. The most difficult part of this expansion is related to the evaluation of the dispersion energy terms. Even though the orientation dependence of the dispersion energy can be derived analytically with relative ease, the coefficients describing the magnitude of the interaction cannot, and considerable effort has been put into calculating the dispersion coefficients as accurately as possible. The intermolecular interaction energy between two homonuclear diatomic molecules with coordinates (R, θ1, θ2, ϕ) can also be represented as a multipole expansion of the type3 V (R ,θ1 ,θ2 ,ϕ) =

■

COMPUTATIONAL DETAILS The quadrupole and hexadecapole moments as well as the dipole−dipole polarizability were calculated using the finite field method,39 which consists of evaluating numerical derivatives of the energy in the presence of a perturbing electric field. The multipole moment tensor components were calculated as the derivative of the energy with respect to the field strength of a quadrupolar or hexadecapolar field. For the evaluation of the polarizability tensor components, second derivatives of the energy are necessary. For numerical stability, however, the evaluation was done by considering perturbations to the dipole moment instead, so that the use of numerical second derivatives could be avoided. For calculation of the numerical derivatives, a seven-point central differences formula was used. The field strengths used in calculating the quadrupole moment, hexadecapole moment, and polarizability tensor components were 1.0 × 10−4 Eh/ea0, 1.0 × 10−6 Eh/ea0, and 1.0 × 10−3 Eh/ea0, respectively (a0 = 0.52917721092 Å, 1 Å = 10−10 m, Eh = 4.35974434 aJ).41 The field strengths were selected by comparing the finite field results with the respective operator eigenvalues at the Hartree−Fock level. The energy evaluations necessary to extract the dispersion energy from the M2−M2 (M = Zn, Cd, Hg) complex were carried out on a grid in (R, θ). The angles θ1 = θ2 = θ were changed from 0° to 90° in steps of 10°. For each angle, the intermolecular distance R was varied from 10 to 35 Å. We have used the experimentally determined bond lengths for the metal dimers: RZn = 4.19 Å, RCd = 4.07 Å, and RHg = 3.63 Å.42 The two-electron integrals and the energies were converged to 16 and 14 decimal places, respectively. The calculated energies converged smoothly without numerical noise. The calculations were performed using the MOLPRO computer program.43 In all calculations, we used the Dunning correlation-consistent polarized valence triple-ζ basis set with small core relativistic pseudopotential (aug-cc-pVTZ-PP) for the metal atoms.44−47 The (n−1)sp electrons of the metal atoms were kept frozen during the calculations. We have used the CCSD(T) level of theory in all calculations unless otherwise stated. The basis set superposition error (BSSE) has been corrected with the counterpoise method.48−50 The counterpoise corrected interaction energy of two species A and B has been calculated as follows:

∑ Vl ,l ,l(R) Gl ,l ,l(θ1,θ2 ,ϕ) 1 2

1 2

(2)

l1, l 2 , l

with the coupled spherical harmonics Gl1,l2,l(θ1,θ2,ϕ) defined as min(l1, l 2)

Gl1, l2 , l(θ1,θ2 ,ϕ) = 4π 2l + 1

∑

⟨l1, m , l 2 , − m|l , 0⟩

m =−min(l1, l 2)

× Ylm1 (θ1,ϕ) Y l−2 m(θ2 ,− ϕ) cos(mϕ)e−2miϕ (3)

where ⟨l1, m, l2, −m|l, 0⟩ is a Clebsch−Gordan coefficient. The coordinates are analogous to those defined in a previous study.29 If R⃗ is a vector connecting the centers of mass of the interacting species, and r1⃗ and r2⃗ are bond vectors of the dimers, then R = |R⃗ |, cos θ1 = r̂1·R̂ , and cos θ2 = r̂2·R̂ . The angle ϕ is the book angle between the planes defined by (r1⃗ , R⃗ ) and (r2⃗ , R⃗ ). There is a close connection between the two expansions presented in eqs 1 and 2. For example, the quadrupole− quadrupole interaction has the same angular dependence as the G2,2,4(θ1,θ2,ϕ) function,28,29 and a similar connection can also be found for the other long-range interaction terms.3 Electrostatic and Induction Energies. The induction and electrostatic energies are simple to model with classical electrodynamics, if the multipole moments and polarizabilities of the interacting species are known.3,38 The electrostatic energy refers to the energy of a molecule in a general nonuniform electric field such as that produced by another molecule some distance R away. However, the electric fields of two interacting molecules can also distort one another, causing an additional interaction called the induction interaction. The geometry dependence of these terms can be calculated by expanding the electric potential in a power series in R. These functions can be calculated once and for all and tabulated for further use. The coefficients describing the magnitude of these interactions depend on the electrodynamic properties of the interacting species and are thus system specific. In this work, we use the multipole moments and polarizabilities computed with the finite field method to calculate these coefficients.39 The formulas describing the interaction can be considerably simplified by taking the symmetry of the system into account.

VAB = EAB − EA(CP) − E B(CP)

(4)

where EAB is the energy of the whole system and E(CP) and A are the energies of the isolated monomers A and B E(CP) B calculated with the AB dimer basis set. 5735



■

Article

RESULTS AND DISCUSSION Multipole Moments and Polarizabilities. To correctly account for the intermolecular dispersion energy, it is important to be able to calculate the electrostatic and induction interactions as accurately as possible, because the energy contributions due to these terms can be comparable in magnitude to the dispersion energy itself. The distance dependence of the leading order electrostatic contribution is EΘΘ ∼ R−5, which is less steep than that of the leading order dispersion term, which decreases as R−6. The other contributions (EΘΦ ∼ R−7, EΘα ∼ R−8, and EΦΦ ∼ R−9) decrease much faster with respect to the distance and are thus less important, although they can still affect the results to a certain degree. As all of these energy terms can be modeled, and the respective multipole moments and polarizabilities efficiently and accurately calculated with the finite field method, we have decided to include them all in the calculations of the dispersion energy. The calculated polarizabilities are in good agreement with other computational results.51,52 It is to be noted that the results for the multipole moments, listed in Table 1, are quite small compared, e.g., with

interesting phenomenon, depicted in Figure 1, is reminiscent of a similar behavior found previously for coinage metal clusters.29 As no numerical fluctuations are evident in the raw energy data, and the electrostatic and induction contributions have been taken explicitly into account, these fluctuations are likely to result from the neglected charge-density anisotropy in the pair potential model used. We have also calculated one test case, the (Zn2)2 cluster at geometry θ = 0°, with the aug-cc-pVQZ-PP basis set. This led to about a 4% increase of the C(ZnZn) 6 coefficient from 208.62 to 217.56 eV Å6. No new values were calculated for the multipole moments or polarizabilities, but changes in these terms would have resulted in only minor corrections, estimated to be less than ±1 eV Å6. We thus believe that a larger basis set would not change the results or conclusions reached in a significant manner. A clear similarity can be noted in the behaviors of the C6 coefficients with respect to the system geometry. The relative magnitudes of the coefficients C(ZnZn) < C(HgHg) < C(CdCd) is, 6 6 6 53 likely due to relativistic effects to polarizability, different from the order of the respective elements in the group (Zn, Cd, Hg). As the shape of an electric field around a molecule can to a satisfactory degree be represented by a multipole expansion, it is reasonable to use this approach to model the variations in the calculated C6 coefficients of Table 2. By least-squares fitting coupled spherical harmonics of eq 3 to model the variations in the calculated C6 coefficients, we are able to accurately account for the geometry dependence of the dispersion coefficients. In addition, the fit naturally divides the interaction into anisotropic (geometry dependent) and isotropic (constant) parts, from which the atomic dispersion coefficients can be obtained. The atomic dispersion coefficients thus obtained are listed in Table 3. Because the elements of group XII are closed shell atoms, we can, to further validate our model, compare the results for the C6 dispersion coefficients with those obtained from a direct calculations on the respective M−M systems. The results in Table 3 show that the dispersion coefficients obtained from the isotropic part of the dispersion energy of our model dimers correspond well with the dispersion coefficients calculated from the atomic system. Orientation Averages of the Dispersion Energy. When dealing with complex systems, one is often more interested in some average of the interactions involved rather than the exact values at a specific geometry. To gain further insight into the properties of the model discussed, we have studied how well can this model be used to model the orientation average of the dispersion interaction. We consider a simple orientation average of the dispersion energy over the angle θ at an intermolecular distance R, E̅disp(R), which is obtained by calculating the average of the dispersion energies Edisp(θ,R) used in our calculations of the dispersion coefficients. The orientation dependence of the leading term of the dispersion interaction between two linear molecules is known and has the form3

Table 1. Calculated Values for the Quadrupole and Hexadecapole Moment and Polarizability Tensor Components Used in This Studya molecule

Θ

Φ

α∥

α⊥

Zn2 Cd2 Hg2

−0.1494 −0.3596 −0.4845

−10.1507 −19.0751 −15.0842

92.26 121.66 91.75

71.37 85.50 63.40

a

All values are given in atomic units. For the quadrupole moment 1 au = 4.48655 × 10−40 C m2, for the hexadecapole moment 1 au = 1.25636 × 10−60 C m4, and for the polarizability 1 au = 1.64878 × 10−41 C2 m2 J−1.41

the results of a previous study regarding coinage metal dimers,29 thus making the electrostatic and induction interactions less important. Dispersion Coefficients. We have calculated the atomic van der Waals C6 dispersion coefficients for the group XII atoms from the respective dimers by using a pair potential model described in previous publications.28,29 It can be seen from the results listed in Table 2 that the calculated dispersion coefficients change as a function of the system geometry. This Table 2. Calculated Atomic C6 Coefficients for Group XII Dimersa θ (deg) 0 10 20 30 40 50 60 70 80 90 av

C(ZnZn) 6 208.62 209.39 207.39 201.05 191.33 181.13 172.91 167.64 164.90 164.04 186.84

(41) (52) (79) (96) (91) (78) (65) (51) (41) (39)

C(CdCd) 6 349.02 350.50 345.12 332.18 313.91 294.83 279.39 269.68 265.11 263.90 306.36

(79) (101) (141) (171) (167) (144) (121) (107) (102) (101)

C(HgHg) 6 280.96 280.08 273.16 261.39 247.33 233.65 222.64 215.70 212.57 211.76 243.92

(146) (135) (126) (118) (102) (84) (70) (59) (53) (51)

(6) Edisp =−

+

γ22 2

C6

1+ R { 6

γ202 2

(3 cos2 θ1 − 1) +

γ022 2

(3 cos2 θ2 − 1)

}

[(2 cos θ1 cos θ2 − sin θ1 sin θ2 cos ϕ)2 − cos2 θ1 − cos2 θ2]

(5)

The coefficients are calculated from different M2−M2 clusters (M = Zn, Cd, Hg). The C6 coefficients are expressed in eV Å6 and the angles in degrees. The uncertainties of the least-squares fits given in parentheses are one-standard errors in the least significant digit. a

where the γ coefficients are, for a given monomer geometry, constants depending on the polarizability anisotropies of the two molecules. The average of this function along any path in 5736



Article

Figure 1. Angle dependence of the C(MM) coefficients for different M2−M2 clusters (M = Zn, Cd, Hg) together with a least-squares fit of the 6 coefficients of the Gl1,l2,l terms defined in eq 3

Table 3. Comparison of Group XII Dispersion Coefficientsa coefficient

isotropic value

average value

atomic

C(ZnZn) 6 C(CdCd) 6 C(HgHg) 6

178.77 290.79 231.35

186.84 306.36 243.92

175.16 283.09 222.18

For additional information, see the Supporting Information. For example, in the leading order this results in ⎡1 r 4 + 4r 2R2 + R4 ⎤ (6) Edisp ⎥ ̅′ = −2C6⎢ 6 + ⎣R (R2 − r 2)5 ⎦

a All values are in eV Å6. The isotropic value refers to the constant term of the fits depicted in Figure 1, and the average value, to the average of the dispersion coefficients reported in Table 2.

It should be noted that this function contains explicitly the bond length r as a parameter, whereas, in eq 5, the dependence of the monomer bond length is implicitly included in the γ constants, because the polarizability depends on the bond length. Another notable difference is that in eq 9 the dependence on the intermolecular distance R is more complicated. In fact, the function in eq 9 has an asymptotic expansion at large R of the form

the coordinate space is thus obviously of the form E̅(6) disp = −c6R−6, where the coefficient c6 is some linear combination of the products of the C6 and the γ coefficients. Higher orders can be treated similarly, and the orientation average becomes −n Edisp ̅ = −∑ cnR

(6)

n

(6) Edisp ̅′

with the cn coefficients being some linear combination of the products of the isotropic Cn coefficients and additional terms due to polarizability anisotropy. On the other hand, considering the pair potential model and a cut of the coordinate space with ϕ = 0, and θ1 = θ2 = θ for two identical linear molecules, the dispersion energy can be represented as

(7)

where r is the monomer bond length. The average value of the term characterized by an integer n ≥ 3 over the range of θ used in this work yields, by integration in the case R > r > 0, ⎡ ⎛ r 2 ⎞⎤ (2n) = −2C2n⎢R−2n + R−2n2F 1⎜n ,n;1; 2 ⎟⎥ Edisp ̅′ ⎝ R ⎠⎦ ⎣ ∞

∑ k=0

Γ(k + n)2 r 2k ⎤ ⎥ (k! )2 Γ(n)2 R2k + 2n ⎥⎦

⎡ ⎛ R2 + r 2 ⎞⎤ = − 2C2n⎢R−2n + (R2 − r 2)−n Pn − 1⎜ 2 ⎟⎥ ⎝ R − r 2 ⎠⎦ ⎣

∞

2

k=0

4R

2 2k ⎤

∑ (k + 1) (2kk++6 2) r

⎥ ⎥⎦ (10)

and thus contains contributions of various R terms as opposed to a clear-cut R−6 distance dependence of E̅ disp. We have compared the ability of the two functions E̅disp and E̅ ′disp to model the calculated average dispersion energy. The comparison was done by carrying out a least-squares fit on the calculated data and comparing the fit residuals and goodness-offit measures for these two models. In both models, only the two leading order terms were considered. Clearly, there is no point in comparing the actual optimized parameter values, because one function represents the total dispersion energy as a function of the intermolecular dispersion coefficient, and the other as a function of the interatomic coefficient. The results of this comparison are shown in Figures 2−5. It can be seen from Figure 2 that E̅d′ isp fits better to the calculated points than E̅ disp. This behavior is due to the higher impact that the close-range points have on the fit of E̅ disp. In fact, the neglect of some of the close-range points yields, as can be seen from Figure 3, a much better fit of E̅ disp. It is noteworthy that the fit of E̅d′ isp remains largely unaffected by this procedure (Figure 4), whereas E̅ disp is much more affected (Figure 5). Overall, these results reveal that, surprisingly, E̅′disp is able to model the actual calculated data more closely at a wider distance range than E̅ disp, although the latter equation should

n

⎡ = −2C2n⎢R−2n + ⎢⎣

⎡1 = −2C6⎢ 6 + ⎢⎣ R

⎞ ⎛4 18r 2 72r 4 = −C6⎜ 6 + 8 + 10 + ...⎟ ⎠ ⎝R R R

′ = −∑ C2n[2R−2n + (r 2 + R2 + 2rR cos θ)−n Edisp + (r 2 + R2 − 2rR cos θ )−n ]

(9)

(8)

where 2F1 is the ordinary hypergeometric function, Γ is the gamma function, and Pn is the Legendre polynomial of order n.54 5737



Article

describe the actual physical phenomenon more accurately. The reasons for this behavior are not entirely clear, but it is possible that the functional form of E̅disp is simply more susceptible to errors in the energy due to numerical factors and the

Figure 5. Least-squares fits of E̅ disp with all points included in the fit (dotted blue curve) and with the first three points neglected in the fit (dashed red curve) to the orientation averaged dispersion energy of the Zn2−Zn2 system.

bond length in E̅ ′disp allow a more concise approximation to some higher order terms not present in E̅disp. It should be noted that, because of the tiny energies involved, this type of fitting can, in general, be ill-behaved and prone to numerical errors, and therefore the nature of the fitting scheme can make a difference between a good fit and a bad one. In addition to the E̅ d′ isp fit described above (model 1), we have also considered a logarithmic fitting scheme (model 2), where ln |Edisp| is fitted for ln R (Figure 6). Under ideal conditions

Figure 2. Least-squares fits of E̅ ′disp (dotted blue curve) and E̅ disp (dashed red curve) to the orientation averaged dispersion energy of the Zn2−Zn2 system.

Figure 3. Least-squares fits of E̅ ′disp (dotted blue curve) and E̅ disp (dashed red curve) to the orientation averaged dispersion energy of the Zn2−Zn2 system with the first three points neglected in the fit. Figure 6. Logarithmic least-squares fit to the orientation averaged dispersion energy for different M2−M2 clusters (M = Zn, Cd, Hg).

these two schemes yield the same results, but in practice the two models differ slightly. To shed light on the general characteristics of these models, we have performed analyses on simulated data of the form of eq 7. Although this simple analysis cannot be Table 4. Dispersion Coefficients Obtained by Orientation Averaginga

Figure 4. Least-squares fits of E̅ ′disp with all points included in the fit (dotted blue curve) and with the first three points neglected in the fit (dashed red curve) to the orientation averaged dispersion energy of the Zn2−Zn2 system.

coefficient

model 1

model 2

atomic

C(ZnZn) 6 C(CdCd) 6 C(HgHg) 6

173.34(119) 279.00(194) 213.12(42)

165.32(26) 266.57(54) 221.97(34)

175.16 283.09 222.18

a

All values are in eV Å6. The numbers in parentheses are one-standard errors in the least significant figure. In the case of model 2, which is nonlinear in the parameters, the errors are calculated on the basis of a linearized approximation, which, however, should provide a reasonable estimate given that the curvature of the parameter space around the fitted values is small in magnitude.

approximations made, e.g., to the electrostatic energy. It can also be the case that if the series representation is truncated, the asymptotic expansion terms and the explicit inclusion of the 5738



Article

(4) Wagner, J. P.; Schreiner, P. R. Nature Utilizes Unusual High London Dispersion Interactions for Compact Membranes Composed of Molecular Ladders. J. Chem. Theory Comput. 2014, 10, 1353−1358. (5) Alonso, J.; Mañanes, A. Long-Range van der Waals Interactions in Density Functional Theory. Theor. Chem. Acc. 2007, 117, 467−472. (6) Cohen, A. J.; Mori-Sánchez, P.; Yang, W. Challenges for Density Functional Theory. Chem. Rev. 2012, 112, 289−320. (7) Burke, K. Perspective on Density Functional Theory. J. Chem. Phys. 2012, 136, 150901. (8) Grimme, S. Density Functional Theory with London Dispersion Corrections. Wiley Interdiscip. Rev.: Comput. Mol. Sci. 2011, 1, 211− 228. (9) Klimeš, J.; Michaelides, A. Perspective: Advances and Challenges in Treating van der Waals Dispersion Forces in Density Functional Theory. J. Chem. Phys. 2012, 137, 120901. (10) Ramalho, J. P. P.; Gomes, J. R. B.; Illas, F. Accounting for van der Waals Interactions Between Adsorbates and Surfaces in Density Functional Theory Based Calculations: Selected Examples. RSC Adv. 2013, 3, 13085−13100. (11) Szalewicz, K. Symmetry-Adapted Perturbation Theory of Intermolecular Forces. Wiley Interdiscip. Rev.: Comput. Mol. Sci. 2012, 2, 254−272. (12) Jansen, G. Symmetry-Adapted Perturbation Theory Based on Density Functional Theory for Noncovalent Interactions. Wiley Interdiscip. Rev.: Comput. Mol. Sci. 2014, 4, 127−144. (13) Korona, T.; Przybytek, M.; Jeziorski, B. Time-Independent Coupled Cluster Theory of the Polarization Propagator. Implementation and Application of the Singles and Doubles Model to Dynamic Polarizabilities and van der Waals Constants. Mol. Phys. 2006, 104, 2303−2316. (14) Norman, P. A Perspective on Nonresonant and Resonant Electronic Response Theory for Time-Dependent Molecular Properties. Phys. Chem. Chem. Phys. 2011, 13, 20519−20535. (15) Mackie, I. D.; DiLabio, G. A. Approximations to Complete Basis Set-Extrapolated, Highly Correlated Non-Covalent Interaction Energies. J. Chem. Phys. 2011, 135, 134318. (16) Hohenstein, E. G.; Sherrill, C. D. Wavefunction Methods for Noncovalent Interactions. Wiley Interdiscip. Rev.: Comput. Mol. Sci. 2012, 2, 304−326. (17) Aradhya, S. V.; Frei, M.; Hybertsen, M. S.; Venkataraman, L. Van der Waals Interactions at Metal/Organic Interfaces at the SingleMolecule Level. Nat. Mater. 2012, 11, 872−876. (18) Chałasinski, G.; Szczesniak, M. M.; Kendall, R. A. Supermolecular Approach to Many-Body Dispersion Interactions in Weak van der Waals Complexes: He, Ne, and Ar Trimers. J. Chem. Phys. 1994, 101, 8860. (19) Elstner, M.; Hobza, P.; Frauenheim, T.; Suhai, S.; Kaxiras, E. Hydrogen Bonding and Stacking Interactions of Nucleic Acid Base Pairs: A Density-Functional-Theory Based Treatment. J. Chem. Phys. 2001, 114, 5149−5155. (20) Du, A. J.; Smith, S. C. Van der Waals-Corrected Density Functional Theory: Benchmarking for Hydrogen-Nanotube and Nanotube-Nanotube Interactions. Nanotechnology 2005, 16, 2118. (21) Williams, R. W.; Malhotra, D. van der Waals Corrections to Density Functional Theory Calculations: Methane, Ethane, Ethylene, Benzene, Formaldehyde, Ammonia, Water, PBE, and CPMD. Chem. Phys. 2006, 327, 54−62. (22) Tekin, A.; Jansen, G. How Accurate is the Density Functional Theory Combined with Symmetry-Adapted Perturbation Theory Approach for CH-π and π-π Interactions? A Comparison to Supermolecular Calculations for the Acetylene-Benzene Dimer. Phys. Chem. Chem. Phys. 2007, 9, 1680−1687. (23) Grimme, S. Semiempirical GGA-type Density Functional Constructed with a Long-Range Dispersion Correction. J. Comput. Chem. 2006, 27, 1787−1799. (24) Grimme, S.; Antony, J.; Ehrlich, S.; Krieg, H. A Consistent and Accurate ab Initio Parametrization of Density Functional Dispersion Correction (DFT-D) for the 94 Elements H-Pu. J. Chem. Phys. 2010, 132, 154104.

considered as a rigorous proof, it suggests some interesting characteristics of the models under study. For example, in general only the coefficients of the lowest order terms can be optimized accurately, and the errors in the coefficients of the higher order terms increase rapidly as the amount of data becomes smaller. Furthermore, the inclusion of higher order terms improves the fitting of the coefficients of lower order terms, with the effect being greatest for “adjacent” terms. Also, the logarithmic model seems to yield better results for the lowest order terms but cannot be termed better overall. The results for the interatomic C6 coefficients, shown in Table 4, are in good comparison with the atomic coefficients.

■

CONCLUSIONS We have modeled dispersion interactions in small clusters of group XII metal atoms with accurate coupled cluster calculations. A pair potential model with a coupled spherical harmonic correction was used to obtain atomic dispersion interaction coefficients from the total intermolecular interaction energy. The lowest order electrostatic and induction terms have been explicitly accounted for, and the necessary multipole moments and polarizabilities have been calculated to a high degree of accuracy with the finite field method. This method of extracting the dispersion coefficients yields results in good agreement with a straightforward atomic calculation, and it is also applicable in the case of open-shell species, where such atomic calculations can be complicated. The coefficients thus obtained could then also be used to describe the dispersion energy in larger systems as a separate contribution to the total interaction energy. We have also investigated the behavior of the orientation averaged dispersion energies. The distance dependence has been modeled by integral averaging of our pair potential model, which results in simple closed form formulas. In the cases studied here, our method results in a more robust fit than a straightforward inverse power law. All of the dispersion coefficients are obtained from ab initio calculations without any experimental parameters.

■

ASSOCIATED CONTENT

S Supporting Information *

The derivation of the formulas describing the orientation averaged dispersion energy in the pair potential model presented in equation eq 8. This material is available free of charge via the Internet at http://pubs.acs.org.

■

AUTHOR INFORMATION

Corresponding Author

*L. Halonen. E-mail: lauri.halonen@helsinki.fi. Phone: +358 2941 50280. Notes

The authors declare no competing financial interest.

■ ■

ACKNOWLEDGMENTS The authors thank the Center for Scientific Computing (CSC) for computational resources. REFERENCES

(1) Kaplan, I. G. Intermolecular Interactions; John Wiley & Sons, Ltd.: New York, 2006. (2) Riley, K. E.; Hobza, P. Noncovalent Interactions in Biochemistry. Wiley Interdiscip. Rev.: Comput. Mol. Sci. 2011, 1, 3−17. (3) Stone, A. J. The Theory of Intermolecular Forces, 2nd ed.; Oxford University Press: Oxford, U.K., 2013. 5739



Article

244. Relativistic Effects in Heavy-Element Chemistry and Physics. In Memoriam Bernd A. Hess (1954−2004). (47) Peterson, K. A.; Puzzarini, C. Systematically Convergent Basis Sets for Transition Metals. II. Pseudopotential-Based Correlation Consistent Basis Sets for the Group 11 (Cu, Ag, Au) and 12 (Zn, Cd, Hg) Elements. Theor. Chem. Acc. 2005, 114, 283−296. (48) Boys, S.; Bernardi, F. The Calculation of Small Molecular Interactions by the Differences of Separate Total Energies. Some Procedures with Reduced Errors. Mol. Phys. 1970, 19, 553−566. (49) van Duijneveldt, F. B.; van Duijneveldt-van de Rijdt, J. G. C. M.; van Lenthe, J. H. State of the Art in Counterpoise Theory. Chem. Rev. 1994, 94, 1873−1885. (50) Tao, F.-M. Bond Functions, Basis Set Superposition Errors and Other Practical Issues with ab Initio Calculations of Intermolecular Potentials. Int. Rev. Phys. Chem. 2001, 20, 617−643. (51) Schautz, F.; Flad, H.-J.; Dolg, M. Quantum Monte Carlo Study of Be2 and Group 12 Dimers M2 (M = Zn, Cd, Hg). Theor. Chem. Acc. 1998, 99, 231−240. (52) Schwerdtfeger, P.; Wesendrup, R.; Moyano, G. E.; Sadlej, A. J.; Greif, J.; Hensel, F. The Potential Energy Curve and Dipole Polarizability Tensor of Mercury Dimer. J. Chem. Phys. 2001, 115, 7401−7412. (53) Seth, M.; Schwerdtfeger, P.; Dolg, M. The Chemistry of the Superheavy Elements. I. Pseudopotentials for 111 and 112 and Relativistic Coupled Cluster Calculations for (112)H+, (112)F2, and (112)F4. J. Chem. Phys. 1997, 106, 3623−3632. (54) Abramowitz, M., Stegun, I. A., Eds. Handbook of Mathematical Functions with Formulas, Graphs, and Mathematical Tables, 10th ed.; Dover Publications: New York, 1972.

(25) Tkatchenko, A.; Scheffler, M. Accurate Molecular van der Waals Interactions from Ground-State Electron Density and Free-Atom Reference Data. Phys. Rev. Lett. 2009, 102, 073005. (26) Zhang, G.-X.; Tkatchenko, A.; Paier, J.; Appel, H.; Scheffler, M. Van der Waals Interactions in Ionic and Semiconductor Solids. Phys. Rev. Lett. 2011, 107, 245501. (27) Tkatchenko, A.; DiStasio, R. A.; Car, R.; Scheffler, M. Accurate and Efficient Method for Many-Body van der Waals Interactions. Phys. Rev. Lett. 2012, 108, 236402. (28) Hänninen, V.; Korpinen, M.; Ren, Q.; Hinde, R.; Halonen, L. An ab Initio Study of van der Waals Potential Energy Parameters for Silver Clusters. J. Phys. Chem. A 2011, 115, 2332−2339. (29) Hatz, R.; Korpinen, M.; Hänninen, V.; Halonen, L. Characterization of the Dispersion Interactions and an ab Initio Study of van der Waals Potential Energy Parameters for Coinage Metal Clusters. J. Phys. Chem. A 2012, 116, 11685−11693. (30) de Heer, W. A. The Physics of Simple Metal Clusters: Experimental Aspects and Simple Models. Rev. Mod. Phys. 1993, 65, 611−676. (31) Flad, H.-J.; Schautz, F.; Wang, Y.; Dolg, M.; Savin, A. On the Bonding of Small Group 12 Clusters. Eur. Phys. J. D 1999, 6, 243−254. (32) Zhao, J. Density-Functional Study of Structures and Electronic Properties of Cd Clusters. Phys. Rev. A: At., Mol., Opt. Phys. 2001, 64, 043204. (33) Moyano, G. E.; Wesendrup, R.; Söhnel, T.; Schwerdtfeger, P. Properties of Small- to Medium-Sized Mercury Clusters from a Combined ab Initio, Density-Functional, and Simulated-Annealing Study. Phys. Rev. Lett. 2002, 89, 103401. (34) Wang, J.; Wang, G.; Zhao, J. Nonmetal-Metal Transition in Znn (n = 2−20) Clusters. Phys. Rev. A: At., Mol., Opt. Phys. 2003, 68, 013201. (35) Gaston, N.; Schwerdtfeger, P. From the van der Waals Dimer to the Solid State of Mercury with Relativistic ab Initio and Density Functional Theory. Phys. Rev. B: Condens. Matter Mater. Phys. 2006, 74, 024105. (36) Buckingham, A. D.; Fowler, P. W.; Hutson, J. M. Theoretical Studies of van der Waals Molecules and Intermolecular Forces. Chem. Rev. 1988, 88, 963−988. (37) Buckingham, A. D. Intermolecular Forces. Pure Appl. Chem. 1970, 24, 123−134. (38) Buckingham, A. D. Advances in Chemical Physics: Intermolecular Forces; John Wiley & Sons, Inc., 1967; Vol. 12, Chapter Permanent and Induced Molecular Moments and Long-Range Intermolecular Forces, pp 107−142. (39) Elking, D. M.; Perera, L.; Duke, R.; Darden, T.; Pedersen, L. G. A Finite Field Method for Calculating Molecular Polarizability Tensors for Arbitrary Multipole Rank. J. Comput. Chem. 2011, 32, 3283−3295. (40) Gelessus, A.; Thiel, W.; Weber, W. Multipoles and Symmetry. J. Chem. Educ. 1995, 72, 505. (41) Mohr, P. J.; Taylor, B. N.; Newell, D. B. CODATA Recommended Values of the Fundamental Physical Constants: 2010. National Institute of Standards and Technology: Gaithesburg, MD, 2010; physics.nist.gov/constants, Retrieved May 15th, 2014. (42) Jules, J. L.; Lombardi, J. R. Transition Metal Dimer Internuclear Distances from Measured Force Constants. J. Phys. Chem. A 2003, 107, 1268−1273. (43) Werner, H.-J.; Knowles, P. J.; Knizia, G.; Manby, F. R.; Schütz, M.; Celani, P.; Korona, T.; Lindh, R.; Mitrushenkov, A.; Rauhut, G.; et al. MOLPRO, Version 2012.1, a Package of ab Initio Programs. 2012; http://www.molpro.net. (44) Dunning, T. H. Gaussian Basis Sets for Use in Correlated Molecular Calculations. I. The Atoms Boron through Neon and Hydrogen. J. Chem. Phys. 1989, 90, 1007. (45) Kendall, R. A.; Dunning, T. H.; Harrison, R. J. Electron Affinities of the First-Row Atoms Revisited. Systematic Basis Sets and Wave Functions. J. Chem. Phys. 1992, 96, 6796. (46) Figgen, D.; Rauhut, G.; Dolg, M.; Stoll, H. Energy-Consistent Pseudopotentials for Group 11 and 12 Atoms: Adjustment to MultiConfiguration Dirac-Hartree-Fock Data. Chem. Phys. 2005, 311, 227− 5740


Dispersion Interactions in Small Zinc, Cadmium, and Mercury Clusters

Recommend Documents