Characterization of the Dispersion Interactions and an ab Initio Study

Oct 26, 2012 - ABSTRACT: Dispersion interactions and van der Waals C6 ... C6 coefficients from a least-squares fit to the interaction energy surface. ...
0 downloads 0 Views 3MB Size
Article pubs.acs.org/JPCA

Characterization of the Dispersion Interactions and an ab Initio Study of van der Waals Potential Energy Parameters for Coinage Metal Clusters Richard Hatz, Markus Korpinen, 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 University of Helsinki, Finland S Supporting Information *

ABSTRACT: Dispersion interactions and van der Waals C6 coefficients are studied in small model systems involving copper, silver, and gold atoms. We investigate a novel method where the intermolecular dispersion interactions are characterized by interatomic C6 coefficients that can be used to formulate the dispersion energy as a separate contribution to the total interaction energy at long distances. We obtain the C6 coefficients from a least-squares fit to the interaction energy surface. Other significant effects to the long-range interaction energy such as multipole−multipole and induction energy are explicitly taken into account. The electronic structure calculations are performed at the coupled cluster level including single, double, and perturbative triple excitations (CCSD(T) level) with triple-ζ basis sets.

1. INTRODUCTION Long-range intermolecular interactions play an important role in a wide range of fields like spectroscopy,1,2 nanoparticles,3−6 molecular crystals,6 and adsorption systems.7 The most important contributions to long-range interactions between molecules are electrostatic, induction, and dispersion interactions. Althoug the former two can be understood in terms of classical electrodynamics, dispersion forces have no real classical analogue but arise because of correlated fluctuations in the charge distributions between interacting atoms or molecules. Although important, the role of dispersion forces is often ignored in calculations of various adsorption systems. These calculations are usually performed using density functional theory (DFT), which lacks the proper description of long-range electron correlation effects.8−10 Even though there have been many attempts to include dispersion interactions into DFT by developing new functionals11−17 or by other schemes,18−22 DFT can still lack some accuracy when dealing with systems where long-range dispersion forces are important. However, the field is currently in a state of rapid change, and new theoretical and computational methods of describing the dispersion interactions are constantly being developed. In addition to the advanced DFT methods, there also exist other theoretical frameworks such as the symmetry-adapted perturbation theory (SAPT)23 and the linear response theory24 as well as many empirical rules,25−28 which could provide the dispersion energy as a distinct contribution to the total interaction energy. Such schemes are, however, often difficult or costly to implement accurately, particularly in open shell species. In this work, we adopt another approach of incorporating dispersion forces into our study of simple van der Waals complexes of coinage metal dimers. This so-called supermolecular method consists of partitioning the total long© 2012 American Chemical Society

range interaction energy into distinct contributions originating from electrostatic, induction, and dispersion forces among other contributions such as many-body interaction terms. In the past, this method has been successfully applied to different kinds of systems in various areas of chemistry, physics, and biology.25,29,30 The basic principle behind the idea, calculating various dispersion coefficients, has received attention in the past,31−33 but the previous studies have almost solely concentrated on the intermolecular interactions between two or several specific molecules thus limiting their applicability in a more general case. Our scheme, on the other hand, aims to construct a more generally applicable model to describe the various long-range interactions between a surface and an adsorbing species. By dividing the dispersion energy into a pairwise additive part and a many-body part, this same approach could be used in different kinds of systems. It would also be a step to the direction of generating a universal model of accounting for long-range electron correlation effects in surface systems. The approach in question has already been implemented by our group in a previous study,34 where promising results were obtained for small silver-containing model systems.

2. LONG-RANGE INTERACTION MODEL 2.1. Intermolecular Interactions. When the interaction between two closed shell molecules is described, the interaction energy can be divided into contributions arising from so-called short-range and long-range effects, which can be further divided Received: July 27, 2012 Revised: October 16, 2012 Published: October 26, 2012 11685

dx.doi.org/10.1021/jp307448n | J. Phys. Chem. A 2012, 116, 11685−11693

The Journal of Physical Chemistry A

Article

into subcategories.31,35 The short-range effects result from the overlap of the wave functions of the interacting molecules and these include the penetration, repulsion, and charge transfer effects caused by the overlap. As the magnitude of the shortrange effects decreases exponentially with respect to the intermolecular distance, these can be ignored if the distance between the interacting molecules is large. At long range, the interaction between the molecules can be conceptually described via the interactions of the charge clouds of the molecules, and the long-range interaction energy Vlong can be further divided into electrostatic, induction, and dispersion contributions Vlong = Velectrostatic + Vinduction + Vdispersion

Figure 1. Jacobi coordinates used in this work. The angles θ1 and θ2 are the angles the dimer bonds form with the line that connects the centers of mass of the two dimers. The angle ϕ12 is the torsional angle that describes the out-of-plane rotation, and the intermolecular distance R is the center-of-mass distance of the dimers. The dimer bond lengths are kept fixed during the calculations.

(1)

The electrostatic and induction energies correspond to the classical long-range interaction terms between charged bodies and the dispersion energy, even though it has no classical counterpart, can be viewed as an attractive contribution to the total energy arising from correlated fluctuations in the electron densities of the interacting species. By using the multipole expansion, the long-range interaction contributions can be further divided into terms that depend on some inverse power of the intermolecular distance R, and the total long-range interaction energy can thus be expressed as Vlong = −∑ CiR−i

Gl1,l2 ,l(θ1 ,θ2 ,ϕ12) = 8π 3/2 min(l1,l 2)

×

× Y l−2 m(θ2 ,−ϕ12) cos(mϕ12)e−2miϕ12

Vinduction = −

1 2

9Θ2 {4(13α⊥ + 17α ) 256R8

+ 3(3α⊥ + 13α )[cos(2θ1) + cos(2θ2)]

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

(4)

where ⟨l1,m,l2,−m|l,0⟩ is a Clebsch−Gordan coefficient. The distance dependence is described by the functions Vl1,l2,l(R), which can be extracted by a least-squares fit to the interaction potential energy surface (abbreviated hereafter as PES). In the fit, we have included all physically meaningful terms of V(R,θ1,θ2,ϕ12) with l1, l2 = 0, 2, 4. 2.3. Induction Interaction. When a molecule is placed in an electric field, electric multipoles may be induced in it because of the distortion of its charge cloud. In the case of interacting molecules, the charge distributions of all of these species are distorted by the electric field and field gradients generated by the other molecules. This interaction between the distorted charge clouds gives rise to induction energy, which can be modeled in terms of classical electrodynamics. The leading order induction energy term in two interacting homonuclear dimers is the quadrupole induced dipole interaction term whose magnitude decreases as R−8. Because of this strong dependence on the intermolecular distance, the induction contributions to the total interaction energy at longrange are small compared to the electrostatic and dispersion contributions. In this work, we have only considered the leading order induction energy calculated as

where the parameters Ci are generally complicated functions depending on the geometry of the system and the multipole moments and polarizabilities of the interacting molecules etc. 2.2. Electrostatic Interaction. Electrostatic interaction is often one of the most important contributions to intermolecular interaction energy at long-range. The general form of the multipole−multipole interaction terms are well-known from classical electrostatics and the magnitude of the interaction can be readily evaluated if the corresponding multipole moments are known.32,36 The leading electrostatic interaction term in this study is the quadrupole−quadrupole interaction, as none of the molecules studied possesses a static dipole moment. We have evaluated the magnitude of the quadrupole and hexadecapole moments of the coinage metal dimers Cu2, Ag2, and Au2 from an analysis of the long-range interactions of the Cu2−H2, Cu2−Cu2, Ag2−H2, and Au2−H2 clusters. The method, which follows that of Diep and Johnson,37,38 has already been described earlier,34 so only a brief account is given here. The interaction between two homonuclear diatomic molecules can be separated into terms describing distance and orientational dependence as V (R ,θ1 ,θ2 ,ϕ12) =

⟨l1 ,m ,l 2 ,−m|l ,0⟩Ylm1 (θ1 ,ϕ12)

∑ m =−min(l1,l 2)

(2)

i

2l + 1 4π

− 10(α⊥ − 3α )[cos(4θ1) + cos(4θ2)]

1 2

(3)

+ 25(α − α⊥)[cos(6θ1) + cos(6θ2)]}

where the coordinates R, θ1, θ2, and ϕ12 describe the orientation of the molecules as depicted in Figure 1. The functions Gl1,l2,l(θ1,θ2,ϕ12) are coupled spherical harmonics describing the angular dependencies in the interaction energy. These are defined as

(5)

between two identical interacting homonuclear dimers, where Θ is the quadrupole moment and α∥ and α⊥ are, respectively, the parallel (αzz) and perpendicular (αxx = αyy) components of the dipole−dipole polarizability tensor. In our calculations, we have used the relativistic values for the polarizability tensor 11686

dx.doi.org/10.1021/jp307448n | J. Phys. Chem. A 2012, 116, 11685−11693

The Journal of Physical Chemistry A

Article

components computed by Saue and Jensen.39 The values used are in Table 1.

hydrogen atom were calculated from the interaction energies of M2−H2 complexes. The C6 coefficients were obtained from a least-squares fit to the dispersion energy, which was obtained from the total interaction energy after subtraction of other significant energy contributions affecting the total interaction energy, such as electrostatic and induction energy terms. The induction energy for the M2−M2 systems was calculated by means of eq 5. The effect of induction energy was found to be small, of the order of a few percent, which is comparable to the inaccuracies originating from the simplicity of the theoretical model used. Thus, we have ignored the induction effects in this study. The basis sets employed in the calculation of the C6 coefficients were the same as those used in the calculation of the multipole moments.

Table 1. Values for the Polarizability Tensor Components Calculated by Saue and Jensen39 Used in This Studya molecule

α∥

α⊥

Cu2 Ag2 Au2

123.8 156.2 114.2

81.9 97.0 65.7

The values are given in atomic units.40 1 au = 1.64878 × 10−41 C2 m2 J−1. a

3. COMPUTATIONAL DETAILS The ab initio calculations necessary to extract the multipole moments were carried out on a grid in (R, θ1, θ2, ϕ12). For the Cu2−H2 and Cu2−Cu2 complexes, we have used 14 and 13 different R values, respectively, ranging from 4.5 to 14 Å (1 Å = 10−10 m). For the Ag2−H2 complex, we have used 11 different values for R from 6 to 12 Å and, for the Au2−H2 complex, we have used 13 different R values from 5 to 12 Å. In all cases, the angles (θ1, θ2, ϕ12) were changed in steps of 30°. The bond lengths of all dimers were kept fixed during the calculations. This was done to simplify the computational treatment and because the changes to the resulting dispersion coefficients due to the relaxation of the bond lengths at such large intermolecular distances were considered negligible. For H2, we have used the vibrationally averaged bond length of the v = 0 state, 0.7668 Å.37 For Cu2, Ag2, and Au2, we have used the experimentally determined bond lengths 2.22,41,42 2.53,43,44 and 2.472 Å,42,45 respectively. The calculations were performed using the MOLPRO computer program.46 In most calculations, we have used the Dunning correlation-consistent polarized valence triple-ζ basis set with small core relativistic pseudopotential47 (aug-cc-pVTZ-PP) for the metal atoms. The only exception was the H2−Cu2 system, where a quadruple-ζ basis set aug-cc-pVQZ-PP was employed for the copper atoms. The basis set for hydrogen atoms was aug-ccpVQZ in all calculations.48−51 The frozen-core approximation has been used for the inner core electrons of the carbon and metal atoms in all 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 in all calculations with the counterpoise method.52−54 The counterpoise corrected interaction energy of two species A and B has been calculated as follows: VAB = EAB − EA(CP) − E B(CP)

4. RESULTS AND DISCUSSION 4.1. Radial Coefficients. At sufficiently long ranges, where short-range effects are negligible, the interaction energy between molecules can be described via long-range interactions such as electrostatic forces along with dispersion and induction interactions. These long-range interactions are described by the radial coefficients Vl1,l2,l(R) each of which describe a different type of contribution to the long-range interaction potential. For example, the V0,0,0 term describes the isotropic component of the interaction, whereas the higher order terms are related to the orientational anisotropy of the long-range interaction, and can be used to extract, e.g., the multipole moments of interest for the systems studied. Figure 2 shows the calculated radial

Figure 2. Calculated radial coefficients for the Cu2−H2 system as a function of the intermolecular distance.

coefficients for the Cu2−H2 system. Complete data of the ab initio calculated radial coefficients for the Cu2−H2, Cu2−Cu2, Ag2−H2, and Au2−H2 clusters together with the corresponding figures are presented in the Supporting Information. We have used the calculated radial coefficients to obtain the quadrupole and hexadecapole moments of the coinage metal dimers. For example, we can evaluate the quadrupole moments of the studied species by noticing that the term G2,2,4(θ1,θ2,ϕ12) calculated from eq 4 has the same angular dependency as the quadrupole−quadrupole interaction energy between two homonuclear dimers.34 This is, of course, no coincidence, but a reflection of the fact that the radial coefficient V2,2,4(R) describes exclusively the quadrupole−quadrupole interaction in the studied systems. Thus, we can evaluate the quadrupole moment from the fit to the PES by using the equation

(6)

E(CP) A

where EAB is the energy of the whole system, and and E(CP) are the energies of the isolated monomers A and B B calculated with the AB dimer basis set. The interatomic van der Waals C6 coefficients were calculated from the counterpoise corrected interaction energies of simple model systems. The two-electron integrals were converged to 16 decimal points and the energies to 14 decimal points. The energies were observed to converge smoothly and no numerical noise was apparent in the energy data. The coefficients C(MM) where M is either Cu, Ag, or Au were 6 calculated from the interaction energies of the corresponding M2−M2 complexes. The coefficients C6(MH) involving the 11687

dx.doi.org/10.1021/jp307448n | J. Phys. Chem. A 2012, 116, 11685−11693

The Journal of Physical Chemistry A

45V2,2,4 4 70

=

Article

reproducing the actual multipole−multipole interaction energies that dominate the total long-range interaction energy. We have used the accurate multipole moment values from our PES surface fit to subtract the lowest order multipole− multipole interaction terms from the total interaction energies of the studied clusters to obtain the van der Waals C6 coefficients for different atom pairs. In the case of the copper dimer, we have used the multipole moment values derived from the Cu2−Cu2 system. 4.2. Interatomic van der Waals C6 Coefficients. 4.2.1. Model Coinage Metal Clusters. Classical considerations are enough to account for electrostatic and induction interaction, but the quantification of dispersion interactions requires special attention. It is known that dispersion effects are approximately pairwise additive meaning that the many-body effects in the dispersion energy are usually small, although recent studies give evidence that this might not be the case for all systems.55 As the different many-body terms are often of different signs, these tend to cancel each other out to a certain degree.32 The induction effects provide another source of nonadditivity to the long-range interaction energy, but their effects can be modeled in classical terms. Furthermore, in our case, the induction effects provide only a minor contribution to the total energy, and as the magnitude of the induction energy depends more strongly on the intermolecular distance than the leading electrostatic or dispersion terms, their contribution is minimal at sufficiently long ranges. We can model the effects of dispersion by calculating the total interaction energy of our small model systems at long intermolecular distances. After subtracting the leading order induction and multipole−multipole interaction energies which in our case consist of the quadrupole−quadrupole, quadrupole−hexadecapole, and the hexadecapole−hexadecapole interaction, we assume that the remaining part of the interaction energy arises predominantly from London type dispersion effects, which we model via the atom pair approximation as

3Θ1Θ2 5

4R

(7)

where Θ1 and Θ2 are the quadrupole moments of the two different molecules in the clusters studied. Similar relations also exist for the computation of higher multipole moments, like hexadecapole moments, which can be calculated from the V2,4,6(R) coefficient. By using the quadrupole moment 0.650 D Å38 (1 D = 3.33564 × 10−30 C m) of H2, we have obtained the results in Table 2 for the multipole moments of the coinage metal dimers. Table 2. Multipole Moments of the Coinage Metal Dimers in Atomic Unitsa molecule

quadrupole moment

hexadecapole moment

Cu2b Cu2c Ag2 Au2

5.72081 5.76314 6.24149 6.4996

56.25985 55.91062 74.55697 117.22882

The atomic unit of the quadrupole moment is 1 au = 4.48655 × 10−40 C m2, and the atomic unit of the hexadecapole moment is 1 au = 1.25636 × 10−60 C m4.40 bCalculated from the Cu2−H2 system. c Calculated from the Cu2−Cu2 system. a

Our values from the potential energy surface fit are in good agreement with values calculated by using a simple finite field method. For the quadrupole moments of Cu2, Ag2, and Au2, we obtained the values 5.8160 ea02, 6.3194 ea02, and 6.5719 ea02, respectively, by using a three-point stencil and a field strength of 1.0 × 10−4 Eh/ea0 (a0 = 0.52918 Å).40 The finite field values for the hexadecapole moments Φ of the coinage metal dimers were ΦCu2 = 56.9952 ea04, ΦAg2 = 74.9375 ea04, and ΦAu2 = 116.4838 ea04. The hexadecapole moments were calculated by using a five-point stencil with a field strength of 1.0 × 10−6 Eh/ ea0. All finite field calculations were done at the CCSD(T)/augcc-pVTZ-PP level. The differences in the multipole moments calculated from the PES fit and by using the finite field method are between 0.5% and 2%, which is a remarkable result especially when we consider the simplicity of the finite field method compared with the PES fit approach and the fact that the method of calculating the multipole moments differs significantly in the two methods. Our approach of calculating the multipole moments from a fit to the potential energy surface is more accurate for our purposes, as the multipole moments are derived from the interaction energies and these are thus capable of better

V = −∑ C6(i)R i−6

The summation is over different atom pairs, and the coefficients C(i) 6 are the interatomic van der Waals dispersion coefficients for the different types of atom pairs, which can be used to model the dispersion effects as a separate contribution to the total interaction energy. Assuming that all atoms lie in the same plane, i.e., ϕ12 = 0, the potential energies of our metal dimer clusters can be written as

⎧ ⎪ 1 V = − C6(MM) × ⎨ ⎪ R2 + 1 R [R R RMM cos(θ1 − θ2) − 2R cos θ2] + θ − 2 cos 1 ⎩ 2 MM MM

{

+

(8)

i

3

}

+

64 2 {[RMM(cos θ1 + cos θ2) + 2R ]2 + RMM (sin θ1 + sin θ2)2 }3

64 1 + 2 1 2 {[RMM(cos θ1 + cos θ2) − 2R ]2 + RMM (sin θ1 + sin θ2)2 }3 R + 2 RMM[RMM − 2R cos θ1 − RMM cos(θ1 − θ2) + 2R cos θ2]

{

⎫ ⎪ ⎬ 3 ⎪ ⎭

}

where M refers to a metal atom and RMM to a metal dimer bond length. The distance R between the dimers and the angles θ1 and θ2 are defined as in Figure 1. Similar relations for other geometries can be derived from eq 8 by means of elementary trigonometry. We have used the method described to calculate the interatomic C6 coefficients C(CuCu) , C(AgAg) , C(AuAu) , C(CuH) , 6 6 6 6 (AgH) (AuH) C6 , and C6 . Table 3 lists the results for the C(CuCu) 6

(9)

coefficient. The results for the other coefficients are in the Supporting Information. From the values listed in Table 3 and corresponding tables listed in the Supporting Information, we can see that the interatomic van der Waals C6 coefficients are somewhat dependent on the geometry of the system used to calculate them. Besides dispersion energy, the most important 11688

dx.doi.org/10.1021/jp307448n | J. Phys. Chem. A 2012, 116, 11685−11693

The Journal of Physical Chemistry A

Article

Table 3. C(CuCu) Coefficients Calculated from Different Cu2− 6 Cu2 Clustersa θ1

θ2

fitting range

CCSD(T)

0 π/18 π/9 π/6 2π/9 5π/18 π/3 7π/18 4π/9 π/2 average

0 π/18 π/9 π/6 2π/9 5π/18 π/3 7π/18 4π/9 π/2

10−35 10−35 10−35 10−35 10−35 10−35 10−35 10−35 10−35 10−35

132.49(28) 139.92(22) 155.30(8) 164.68(10) 158.97(19) 142.31(20) 128.57(22) 128.23(37) 137.65(64) 143.35(77) 143.15

the interaction energies at long intermolecular distances where the absolute value of the energy is small, the computed errors introduced by the approximations done in the electronic structure calculations could affect the accuracy of the fit and thus the calculated coefficients. For example, the full computational inclusion of many-body effects in the dispersion forces would require a high-order coupled cluster method. Although the exact computational and theoretical inclusion of many-body phenomena in the calculation of our interatomic C6 coefficients would pose an interesting research topic for the future, it is questionable whether the much increased computational costs would be warranted considering the possible gains and intended use of our method as a quick and simple way to gauge the effects of dispersion interactions in large systems. Despite the fact that high-order many-body terms were not explicitly included in our theoretical description, we have been able to model the behavior of the interatomic C6 coefficients accurately as a function of the geometry of the system, as can be seen from Figures 3−5. The differences were accounted for by

The C(CuCu) coefficients are expressed in eV Å6, the fitting ranges in 6 Å, and the angles in radians. The uncertainties given in parentheses are one-standard errors in the least significant digit. a

contribution to the long-range interaction energy is the quadrupole−quadrupole energy, and the accurate subtraction of this electrostatic interaction term is vital to the accurate calculation of the interatomic C6 coefficients. We have been able to rule out possible inaccuracies in the calculated quadrupole moments as sources of error in the van der Waals coefficients by calculating the coefficients from dimer orientations where the quadrupole−quadrupole interaction contribution vanishes. For example, it can be seen from Table 5 of the Supporting Information that the differences in the C(AgAg) persist even in cases when the orientations of the 6 silver dimers are chosen such that the effect of quadrupole− quadrupole interaction is nonexistent. We have also considered the possibility of uncertainties in the higher order electrostatic or induction interaction terms affecting the accuracy of the calculated C6 coefficients. Although we observed that the inclusion of higher order terms helped to mitigate the differences in the calculated C6 terms to some degree, the effect of possible uncertainties in these terms is not large enough to account for the observed discrepancies. Because of the simplicity of our model and the fact that the uncertainties in the most significant long-range interaction terms are unable to account for the differences in the C6 coefficients, it seems plausible that these differences are related to the neglected higher order many-body terms in the theoretical description of the dispersion forces. We have investigated the possibility of accounting for the leading order many-body interaction energy term explicitly by including Axilord−Teller−Muto56,57 three-body interaction terms (ATM terms) into the theoretical model used, but it was found that these terms were unable to account for the differences observed. The results are not surprising because there have been many other studies indicating that that ATM terms, although useful in many situations,58−66 incorrectly account for the many-body interactions in some systems.63,67−69 Moreover, the perturbation theory approach used to derive the ATM terms is strictly valid only when the distance between all interacting atoms is largea condition that is not satisfied in our simple model systems, although the theory has also been used regardless of interatomic distance because of its apparent success and convenience.70 Besides the many-body interaction terms, other factors can also lead to the discrepancies observed in the dispersion coefficients. As the C6 coefficients are calculated from a fit to

Figure 3. C(CuCu) coefficients for different Cu2−Cu2 clusters. The 6 intermolecular distance is fixed to R = 10 Å. The curve depicts a fit of the form shown in eq 3 to the C6 coefficients.

Figure 4. C(AgAg) coefficients for different Ag2−Ag2 clusters. The 6 intermolecular distance is fixed to R = 10 Å. The curve depicts a fit of the form shown in eq 3 to the C6 coefficients.

fitting the dispersion coefficients to a function of the form of eq 3. Explicitly, the function used can be written as (MM) V disp (θ1 ,θ2 ,ϕ12 ,R ) =

∑ Cl(MM) ,l ,l (R ) G l ,l ,l (θ1 ,θ2 ,ϕ12) 1 2

1 2

l1,l 2 ,l

(10)

where the coefficients C(MM) l2,l2,l (R) are the fitting parameters and M is a metal atom. In our least-squares fit, we have included the 11689

dx.doi.org/10.1021/jp307448n | J. Phys. Chem. A 2012, 116, 11685−11693

The Journal of Physical Chemistry A

Article

we have, for example, shown that this method could be used to gain insight on the adsorption of a silver cluster on a graphene surface. However, it should also be noted that recent studies indicate that some hybrid inorganic−organic systems might require a more in-depth theoretical and computational treatment.71 4.2.2. Toward Cluster−Surface Interactions: Silver Dimer− Naphthalene Model. This study of interatomic van der Waals dispersion coefficients was partly motivated by our interest in the dispersion forces between a metal cluster and a graphene surface. In a previous paper,34 we noticed by examining the silver dimer−benzene complex that it might be possible to model the dispersion interaction between a metal cluster and a surface by using our quick and simple method of interactomic C6 coefficients. We have also calculated the C(AgC) coefficient 6 from a larger silver dimer−naphthalene complex to compare the results with our previous calculations with the benzene complex. We have used a previously optimized geometry72 at the B3LYP/6-311++G(d,p) level for the naphthalene molecule described in Figure 6 and three different models (Figures 7−9)

Figure 5. C(AuAu) coefficients for different Au2−Au2 clusters. The 6 intermolecular distance is fixed to R = 10 Å. The curve depicts a fit of the form shown in eq 3 to the C6 coefficients.

leading order terms allowed by symmetry, which in our model systems amount to the terms with indices (l1, l2, l) = (0, 0, 0), (2, 0, 2), (2, 2, 4), (2, 4, 6), and (4, 4, 8). The most important contribution for the study of interatomic C6 coefficients is the C(MM) 0,0,0 (R) term, which describes the isotropic part of the interaction and is the only term that is present between two interacting spherically symmetric atoms, and can thus be related to the interatomic C6 coefficients. The other terms in the fit could be interpreted as describing the residual anisotropies caused by the use of the simplistic pair potential model. The fact that the interatomic C6 coefficients derived from the fit, reported in Table 4, coincide well with the average values Table 4. Interatomic C6 Coefficients in eV Å6 Derived from the Isotropic Part of the V(MM) disp Function and Comparison with the Average Values from Table 3 and Tables 5 and 6 of the Supporting Informationa coefficient

value from fit

average value

C(CuCu) 6 C(AgAg) 6 C(AuAu) 6

141.72(4) 231.70(17) 215.53(33)

143.15 234.67 216.96

Figure 6. Geometry of the naphthalene molecule used in the calculations.

a

The uncertainties given in parentheses are one-standard errors in the least significant digit.

reported in Table 3 and Tables 5 and 6 of the Supporting Information is reassuring, because it shows that even though there are geometry dependent discrepancies in the interatomic C6 coefficients, the orientational average or the isotropic part can still be used to estimate the dispersion energy in large systems. The slight differences in the values of the van der Waals coefficients do not much affect the applicability of our method in the computational quantitative determination of dispersion forces as a separate contribution to the interaction energy. The strength of our method lies in the point that it can be used to model the total dispersion energy simply in large systems, where the small geometrical dependencies tend to cancel each other out so that the total dispersion energy can be gauged by including the many-body effects in an average fashion. Thus, this method could be used to estimate the dispersion contribution to the total interaction energy in systems like molecular clusters, biomolecules, or in reactions of a metal cluster with a solid surface, where the effect of dispersion forces is important yet difficult to calculate accurately with the methods currently used to study them. In a previous paper,34

Figure 7. Silver−naphthalene complex 1.

coefficient as shown in Table 5. to calculate the average C(AgC) 6 To compare the average C(AgC) coefficients calculated from the 6 benzene and naphthalene complexes, we have, in the case of the naphthalene complexes, used the same electronic structure calculation methods (Møller−Plesset perturbation theory, MP2 and its spin component scaled version, SCS-MP273) and the silver dimer quadrupole moment (ΘAg2 = 8.132 51 D Å) as in the previous benzene complex calculations. The MP2 values are given mainly for comparison because the MP2 method is known to overestimate the interaction energy of complex systems74−76 leading to too large C6 coefficientsa problem 11690

dx.doi.org/10.1021/jp307448n | J. Phys. Chem. A 2012, 116, 11685−11693

The Journal of Physical Chemistry A

Article

Even though the dispersion coefficients calculated from the different models differ somewhat from each other, it is interesting to note that the average coefficient C(AgC) = 77.75 6 eV Å6 is within 10% from the average coefficient C(AgC) = 71.10 6 eV Å6 calculated previously from the benzene model complexes. The exhaustive investigation of the interactomic C(AgC) coefficient as a function of system size and at the bulk 6 limit is outside the scope of the present study, and requires further research. Nevertheless, the fact that we have been able to obtain similar results for the C(AgC) coefficient from both 6 benzene and naphthalene complexes is a further indication that our model could be used to estimate the magnitude of dispersion interactions in large systems as a separate contribution18,25,29,30,81 to the total interaction energy. For example, the C(AgC) coefficient is especially interesting, because 6 the dispersion forces have been shown to be important in the computational modeling of metal clusters adsorbed to a graphene surface.82

Figure 8. Silver−naphthalene complex 2.

5. CONCLUSIONS The role of van der Waals dispersion forces in large systems has received an ever increasing amount of attention during the past few years. Dispersion effects are, however, still occasionally overlooked in computational studies of, e.g., adsorption systems, molecular complexes, metal clusters, and biomolecules although these forces can significantly contribute to the energetics of the system. In fact, weak interactions are often ill characterized in most common DFT functionals usually employed in the study of large systems. Even though the field has seen some advancement in the form of new functionals and other, rather involved, means to account for long-range interactions within the DFT formalism, a simple generally applicable scheme of accounting for dispersion interactions in large systems is called for. The obvious solution of using a more sophisticated level of theory to describe the system is unfortunately often inapplicable, though, because of the computational resources available. Our method of quantifying the dispersion interactions as a separate contribution to the total interaction energy, on the other hand, is simple and applicable to many different kinds of systems once the interatomic van der Waals C6 coefficients are known. Thus, the problem of having to account for the dispersion interactions implicitly via complicated schemes can therefore be reduced to calculate the dispersion interactions explicitly via the interatomic dispersion coefficients. We have calculated dispersion coefficients between several different pairs of atoms involving the coinage metals Cu, Ag, and Au, which provide an interesting theoretical test case considering the exceptional amount of scientific interest directed toward the coinage metals and metal clusters in general over the past decades.83−86 We have calculated the van der Waals coefficients on the basis of rigorous ab initio methods without the need for experimental parameters such as the Slater−Kirkwood effective number87 among other parameters often used to assess the significance of dispersion interactions in complex systems.

Figure 9. Silver−naphthalene complex 3.

Table 5. Average C(AgC) Coefficient Calculated from the 6 Different Silver Dimer−Naphthalene Models in eV Å6 a model model 1 model 2 model 3 average C(AgC) 6

fitting range

MP2

SCS-MP2

8−15 8−20 8−15 8−20 8−15 8−20 8−15 8−20

111.08(176)

85.76(143) 85.77(166) 95.91(149) 95.92(172) 51.56(101) 51.55(116) 77.74 77.75

121.22(181) 74.52(125) 102.27

a The fitting ranges are in Å. The uncertainties given in parentheses are one-standard errors in the least significant digit.

which the SCS-MP2 method largely corrects.73,77 Our method of accurately calculating the quadrupole moment from the fitted radial coefficients cannot be used for naphthalene, but fortunately it has been shown that the Hartree−Fock (HF) method produces accurate quadrupole moments for benzene and naphthalene.78,79 We have used the values Θzz = −9.2015 D Å, Θxx = 4.4947 D Å, and Θyy = 4.7068 D Å for the quadrupole moment tensor elements of the naphthalene molecule calculated with the Gaussian09 program80 at the HF/aug-cc-pVQZ level.



ASSOCIATED CONTENT

S Supporting Information *

Additional radial coefficients, including figures of radial coefficients as a function of the intermolecular distance, and intermolecular dispersion coefficients not given in the article. 11691

dx.doi.org/10.1021/jp307448n | J. Phys. Chem. A 2012, 116, 11685−11693

The Journal of Physical Chemistry A

Article

(40) Codata recommended values of the fundamental physical constants: 2010. Mohr, P. J.; Taylor, B. N.; Newell, D. B. National Institute of Standards and Technology, Gaithesburg, MD. (41) Page, R. H.; Gudeman, C. S. J. Chem. Phys. 1991, 94, 39. (42) Haynes, W. M., Ed.; CRC Press/Taylor and Francis, Boca Raton, FL, 2012; Chapter Spectroscopic Constants of Diatomic Molecules. (43) Ran, Q.; Schmude, W.; Gingerich, K. A.; Wilhite, D. W. J. Phys. Chem. 1993, 97, 8535. (44) Krämer, H.; Beutel, V.; Weyers, K.; Demtröder, W. Chem. Phys. Lett. 1992, 193, 331. (45) Simard, B.; Hackett, P. A. J. Mol. Spectrosc. 1990, 412, 310. (46) Werner, H.-J.; Knowles, P. J.; Manby, F. R.; Schütz, M.; Celani, P.; Knizia, G.; Korona, T.; Lindh, R.; Mitrushenkov, A.; Rauhut, G.; Adler, T. B.; Amos, R. D.; Bernhardsson, A.; Berning, A.; Cooper, D. L.; Deegan, M. J. O.; Dobbyn, A. J.; Eckert, F.; Goll, E.; Hampel, C.; Hesselmann, A.; Hetzer, G.; Hrenar, T.; Jansen, G.; Köppl, C.; Liu, Y.; Lloyd, A. W.; Mata, R. A.; May, A. J.; McNicholas, S. J.; Meyer, W.; Mura, M. E.; Nicklass, A.; Palmieri, P.; Pflüger, K.; Pitzer, R.; Reiher, M.; Shiozaki, T.; Stoll, H.; Stone, A. J.; Tarroni, R.; Thorsteinsson, T.; Wang, M.; Wolf, A. Molpro, version 2010.1, a package of ab initio programs, 2010. (47) Figgen, D.; Rauhut, G.; Dolg, M.; Stoll, H. Chem. Phys. 2005, 311, 227. (48) Dunning, T. H. J. Chem. Phys. 1989, 90, 1007. (49) Kendall, R. A.; Dunning, T. H.; Harrison, R. J. J. Chem. Phys. 1992, 96, 6796. (50) Feller, D. J. Comput. Chem. 1996, 17, 1571. (51) Schuchardt, K. L.; Didier, B. T.; Elsethagen, T.; Sun, L.; Gurumoorthi, V.; Chase, J.; Li, J.; Windus, T. L. J. Chem. Inf. Model. 2007, 47, 1045. (52) Boys, S. F.; Bernardi, F. Mol. Phys. 1970, 19, 553. (53) van Duijneveldt, F. B.; van Duijneveldt-van de Rijdt, J. G. C. M.; van Lenthe, J. H. Chem. Rev. 1994, 94, 1873. (54) Tao, F. Int. Rev. Phys. Chem. 2001, 20, 617. (55) von Lilienfeld, O. A.; Tkatchenko, A. J. Chem. Phys. 2010, 132, 234109. (56) Axilord, P. M.; Teller, E. J. Chem. Phys. 1943, 11, 299. (57) Muto, Y. Proc. Phys.-Math. Soc. Jpn. 1943, 17, 629. (58) Dymond, J. H.; Rigby, M.; Smith, E. B. J. Chem. Phys. 1965, 42, 2801. (59) Barker, J. A.; Pompe, A. Aust. J. Chem. 1968, 21, 1683. (60) Doran, M. B.; Zucker, I. J. J. Phys. C 1971, 4, 307. (61) Barker, J. A.; Watts, R. O.; Lee, J. K.; Schafer, T. P.; Lee, Y. T. J. Chem. Phys. 1974, 61, 3081. (62) Etters, R. D.; Danilowicz, R. J. Chem. Phys. 1979, 71, 4767. (63) Elrod, M. J.; Saykally, R. J. Chem. Rev. 1994, 94, 1975. (64) Anta, J. A.; Lomba, E. Phys. Rev. E 1997, 55, 2707. (65) Bukowski, R.; Szalewicz, K. J. Chem. Phys. 2000, 114, 9518. (66) Leonhard, K.; U. K. Deiters, U. K. Mol. Phys. 2000, 98, 1603. (67) Egelstaff, P. A.; Teitsma, A. Phys. Rev. A 1980, 22, 1702. (68) Kalos, M. H.; Lee, M. A.; Wihitlock, P. A.; Chester, G. V. Phys. Rev. B 1981, 24, 115. (69) Donchev, A. G. J. Chem. Phys. 2006, 125, 074713. (70) Maitland, G. C.; Rigby, M.; Smith, E. B.; Wakeham, W. A. Intermolecular Forces: Their Origin and Determination; Clarendon Press: Oxford, U.K., 1981. (71) Ruiz, V. G.; Liu, W.; Zojer, E.; Scheffler, M.; Tkatchenko, A. Phys. Rev. Lett. 2012, 108, 146103. (72) Kapinus, V. A. Photophysical Properties of Protonated Aromatic Hydrocarbons Ph.D. thesis, California Institute of Technology, 2005. (73) Grimme, S. J. Chem. Phys. 2003, 118, 9095. (74) Hobza, P.; Selzle, H. L.; Schlag, E. W. J. Phys. Chem. 1996, 18790. (75) Sinnokrot, M. O.; Valeev, E. F.; Sherrill, C. D. J. Am. Chem. Soc. 2002, 124, 10887. (76) Tsuzuki, S.; Kazumasa, H.; Uchimaru, T.; Mikami, M. J. Chem. Phys. 2004, 120, 647. (77) Piacenza, M.; Grimme, S. ChemPhysChem 2005, 6, 1554.

This information is available free of charge via the Internet at http://pubs.acs.org.



AUTHOR INFORMATION

Notes

The authors declare no competing financial interest.



REFERENCES

(1) Kouzov, A. P.; Chrysos, M.; Rachet, F.; Egorova, N. I. Phys. Rev. A 2006, 74, 012723. (2) Chrysos, M.; Rachet, F.; Egorova, N. I.; Kouzov, A. P. Phys. Rev. A 2007, 75, 012707. (3) Pacheco, J. M.; Prates Ramalho, J. P. Phys. Rev. Lett. 1997, 79, 3873. (4) Grifalco, L. A.; Hodak, M. Phys. Rev. B 2002, 65, 125404. (5) Amovilli, C.; March, N. H. Carbon 2005, 43, 1634. (6) Gibbs, G. V.; Crawford, T. D.; Wallace, A. F.; Cox, D. F.; Parrish, R. M.; Hohenstein, E. G.; Sherrill, C. D. J. Phys. Chem. A 2011, 115, 12933. (7) Dappe, Y. J.; Ortega, J.; Flores, F. Weak Chemical Interaction and van der Waals Forces: A Combined Density Functional and Intermolecular Perturbation Theory - Application to Graphite and Graphitic Systems. Lecture Notes in Physics; Springer: Berlin/ Heidelberg, 2010; Vol. 795, p 45. (8) Alonso, J. A.; Mañanes, A. Theor. Chem. Acc. 2007, 117, 467. (9) Cohen, A. J.; Mori-Sánchez, P.; Yang, W. Chem. Rev. 2012, 112, 289. (10) Burke, K. J. Chem. Phys. 2012, 136, 150901. (11) Andersson, Y.; Langreth, D. C.; Lunqvist, B. I. Phys. Rev. Lett. 1996, 76, 102. (12) Dion, M.; Rydberg, H.; Schröder, E.; Langreth, D. C.; Lundqvist, B. Phys. Rev. Lett. 2004, 92, 246401. (13) Becke, A. D.; Johnson, E. R. J. Chem. Phys. 2005, 122, 154104. (14) Becke, A. D.; Johnson, E. R. J. Chem. Phys. 2005, 123, 154101. (15) Johnson, E. R.; Becke, A. D. J. Chem. Phys. 2005, 123, 024101. (16) Becke, A. D.; Johnson, E. R. J. Chem. Phys. 2006, 124, 014104. (17) Johnson, E. R.; Becke, A. D. J. Chem. Phys. 2006, 124, 174104. (18) Wu, Q.; Yang, W. J. Chem. Phys. 2002, 116, 515. (19) Grimme, S. J. Comput. Chem. 2006, 27, 1787. (20) Tkatchenko, A.; Scheffler, M. Phys. Rev. Lett. 2009, 102, 073005. (21) Zhang, G.-X.; Tkatchenko, A.; Paier, J.; Appel, H.; Scheffler, M. Phys. Rev. Lett. 2011, 107, 245501. (22) Tkatchenko, A.; DiStasio, R. A., Jr.; Car, R.; Scheffler, M. Phys. Rev. Lett. 2012, 108, 236402. (23) Jeziorski, B.; Moszynski, R.; Szalewicz, K. Chem. Rev. 1994, 94, 1887. (24) Norman, P. Phys. Chem. Chem. Phys. 2011, 13, 20519. (25) Elstner, M.; Hobza, P.; Frauenheim, T.; Suhai, S.; Kaxiras, E. J. Chem. Phys. 2001, 114, 5149. (26) Amar, J. Phys. Rev. B 2003, 67, 165425. (27) Surong, Q.; Zhao, Y.; Jing, X.; Li, X.; Su, W. Int. J. Quantum Chem. 2004, 100, 293. (28) Grimme, S. J. Chem. Phys. 2006, 124, 034108. (29) Du, A. J.; Smith, S. C. Nanotechnology 2005, 16, 2118. (30) Williams, R. W.; Malhotra, D. Chem. Phys. 2006, 327, 54. (31) Buckingham, A. D.; Fowler, P. W.; Hutson, J. M. Chem. Rev. 1988, 88, 963. (32) Stone, A. J. The Theory of Intermolecular Forces; Clarendon Press: Oxford, U.K., 1996. (33) Tkatchenko, A.; Romaner, L.; Hofmann, O. T.; Zojer, E.; Ambrosch-Draxl, C.; M., S. MRS Bull. 2010, 35, 435. (34) Hänninen, V.; Korpinen, M.; Ren, Q.; Hinde, R.; Halonen, L. J. Phys. Chem. A 2011, 115, 2332. (35) Buckingham, A. D. Pure Appl. Chem. 1970, 24, 123. (36) Buckingham, A. D. Adv. Chem. Phys. 1967, 12, 107. (37) Diep, P.; Johnson, J. K. J. Chem. Phys. 2000, 112, 4465. (38) Diep, P.; Johnson, J. K. J. Chem. Phys. 2000, 113, 3480. (39) Saue, T.; Jensen, H. J. A. J. Chem. Phys. 2003, 118, 522. 11692

dx.doi.org/10.1021/jp307448n | J. Phys. Chem. A 2012, 116, 11685−11693

The Journal of Physical Chemistry A

Article

(78) Heard, G. L.; Boyd, R. J. Chem. Phys. Lett. 1997, 277, 252. (79) Hernández-Trujillo, J.; Vela, A. J. Phys. Chem. 1996, 100, 6524. (80) Frisch, M. J.; Trucks, G. W.; Schlegel, H. B.; Scuseria, G. E.; Robb, M. A.; Cheeseman, J. R.; Scalmani, G.; Barone, V.; Mennucci, B.; Petersson, G. A.; Nakatsuji, H.; Caricato, M.; Li, X.; Hratchian, H. P.; Izmaylov, A. F.; Bloino, J.; Zheng, G.; Sonnenberg, J. L.; Hada, M.; Ehara, M.; Toyota, K.; Fukuda, R.; Hasegawa, J.; Ishida, M.; Nakajima, T.; Honda, Y.; Kitao, O.; Nakai, H.; Vreven, T.; Montgomery, J. A., Jr.; Peralta, J. E.; Ogliaro, F.; Bearpark, M.; Heyd, J. J.; Brothers, E.; Kudin, K. N.; Staroverov, V. N.; Kobayashi, R.; Normand, J.; Raghavachari, K.; Rendell, A.; Burant, J. C.; Iyengar, S. S.; Tomasi, J.; Cossi, M.; Rega, N.; Millam, J. M.; Klene, M.; Knox, J. E.; Cross, J. B.; Bakken, V.; Adamo, C.; Jaramillo, J.; Gomperts, R.; Stratmann, R. E.; Yazyev, O.; Austin, A. J.; Cammi, R.; Pomelli, C.; Ochterski, J. W.; Martin, R. L.; Morokuma, K.; Zakrzewski, V. G.; Voth, G. A.; Salvador, P.; Dannenberg, J. J.; Dapprich, S.; Daniels, A. D.; Farkas, Ö .; Foresman, J. B.; Ortiz, J. V.; Cioslowski, J.; Fox, D. J. Gaussian09, Revision A.1; Gaussian Inc.: Wallingford, CT, 2009. (81) Grimme, S. J. Comput. Chem. 2004, 25, 1463. (82) Jalkanen, J.-P.; Halonen, M.; Fernández-Torre, D.; Laasonen, K.; Halonen, L. J. Phys. Chem. A 2007, 111, 12317. (83) Pyykkö, P. Angew. Chem., Int. Ed 2004, 43, 4412. (84) Pyykkö, P. Inorg. Chim. Acta 2005, 358, 4113. (85) Pyykkö, P. Chem. Soc. Rev. 2008, 37, 1967. (86) Lu, Y.; Chen, W. Chem. Soc. Rev. 2012, 41, 3594. (87) Slater, J. C.; Kirkwood, J. G. Phys. Rev. 1931, 37, 682.

11693

dx.doi.org/10.1021/jp307448n | J. Phys. Chem. A 2012, 116, 11685−11693