Basis Set and Method Dependence in Quantum Theory of Atoms in

Nov 4, 2010 - Department of Quantum Chemistry, Nicolaus Copernicus UniVersity, 7-Gagarina St., PL-87 100 Torun, Poland, and. Department of Theoretical...
1 downloads 0 Views 3MB Size
Download PDF
12498

J. Phys. Chem. A 2010, 114, 12498–12505

Basis Set and Method Dependence in Quantum Theory of Atoms in Molecules Calculations for Covalent Bonds Mirosław Jabłon´ski*,† and Marcin Palusiak‡ Department of Quantum Chemistry, Nicolaus Copernicus UniVersity, 7-Gagarina St., PL-87 100 Torun´, Poland, and Department of Theoretical and Structural Chemistry, UniVersity of Ło´dz´, 12-Tamka St., PL-91 403 Ło´dz´, Poland ReceiVed: July 20, 2010; ReVised Manuscript ReceiVed: October 8, 2010

The influence of various small- and medium-size basis sets used in Hartree-Fock (HF) and density functional theory (DFT)/B3LYP calculations on results of quantum theory of atoms in molecules based (QTAIM-based) analysis of bond parameters is investigated for several single, double, and triple covalent bonds. It is shown that, in general, HF and DFT/B3LYP methods give very similar QTAIM results with respect to the basis set. The smallest 6-31G basis set and DZ-quality basis sets of Dunning type lead to poor results in comparison to those obtained by the most reliable aug-cc-pVTZ. On the contrary, 6-311++G(2df,2pd) and in a somewhat lesser extent 6-311++G(3df,3pd) basis sets give satisfactory values of QTAIM parameters. It is also demonstrated that QTAIM calculations may be sensitive for the method and basis set in the case of multiple and more polarized bonds. Introduction As opposed to the wave function quantum mechanics, the quantum theory of atoms in molecules (QTAIM) introduced by Bader1-6 gives the opportunity to have an insight into a region of a system. This opportunity meets the interest of most of chemists who wish to have a theoretical tool to study a small part of a molecule only, instead of dealing with the total energy of a whole system. Nowadays QTAIM is a very powerful method in investigating inter- and intramolecular interactions.7-13 In the QTAIM approach a system is divided to subsystems by means of a physically well-defined surface which is established in such a way that there must not be a flux in the electron b F) through the surface. Followdensity gradient vector field (∇ ing the gradient vector of F from some initial point one b F ) f0. Of terminates at a critical point (CP), a point where ∇ particular interest is a critical point in which the electron density is a minimum with respect to the direction of a line connecting two nuclei and a maximum with respect to all directions perpendicular to this line, thus being a saddle point with the (3,-1) curvature. A critical point with such characteristics is called a bond critical point (BCP). BCPs are of particular interest since they can be used in the recognition of chemical bonds between atoms and then in studies of the character of these interatomic interactions. Most of QTAIM calculations are performed with the use of rather small basis sets. Therefore, one may ask about the dependence of values of QTAIM parameters on both the basis set and the method used in QTAIM calculations, i.e., about the stability of values of QTAIM parameters with respect to the change of both the basis set and the method used in QTAIM calculations. The present paper is a continuation of our recent studies12 of this issue. Previously, we were investigating hydrogen bonds of the CH · · · N type.12 In this paper we * To whom correspondence should be addressed. E-mail: teojab@ chem.uni.torun.pl. Tel.: +48 (56) 6114695. Fax: +48 (56) 6542477. † Nicolaus Copernicus University. ‡ University of Ło´dz´.

investigate several covalent bonds, which can be grouped into either single, double, or triple or either nonpolarized or polarized. Methodology To study the stability of QTAIM calculations with respect to the basis set and the method, we chose several small molecules possessing single, double, and triple covalent bonds: (i) ethane (H3CsCH3), fluoromethane (H3CsF), and methanol (H3CsOH), (ii) ethene (H2CdCH2) and formaldehyde (H2CdO), (iii) acetylene (HCtCH), nitrogen molecule (NtN), and hydrogen cyanide (HCtN), where in parentheses we indicated the bond of interest, i.e., that one for which BCP’s values of QTAIM parameters were computed. Geometry optimizations of studied molecules have been performed at the B3LYP/aug-cc-pVQZ level of theory. However, it should be mentioned that in fact any other method could have been used for the geometry optimization since the aim of these studies was to compare values of QTAIM parameters obtained by means of different approximations for a system with the same geometry, whatever the source of the geometry is. Various QTAIM parameters which are commonly used in the description of chemical bonding were then computed using a large set of basis sets of Pople14 and Dunning15,16 type. Smallsize basis sets of Pople type have also been used in these studies since their use may be the only option in the case of much more time-consuming calculations, e.g., for very bulky molecular systems. Results obtained by means of Pople-type basis sets (or the smallest DZ-quality of Dunning type) will also be compared to those obtained by means of more reliable12 ccpVTZ and particularly aug-cc-pVTZ basis sets. In the present calculations we do not use Dunning-type basis sets with the cardinal number larger than 3 (TZ) since the wave function file (wfn) is not generated if the angular momentum quantum number l is greater than 3; i.e., there is a technical limitation for the usage of f-type Gaussian-type orbitals. In our recent studies of QTAIM stability for the CH · · · N hydrogen bond12 we have also used d-aug-cc-pVTZ as well as cc-pVQZ and augcc-pVQZ basis sets with g-functions removed in the case of

10.1021/jp106740e  2010 American Chemical Society Published on Web 11/04/2010

QTAIM Calculations for Covalent Bonds QZ quality basis sets. However, those calculations have shown that aug-cc-pVTZ basis set leads to values of QTAIM parameters close to those obtained by means of larger Dunning-type basis sets. Thus, in the present calculations for relevant covalent bonds, values obtained by means of aug-cc-pVTZ basis set are assumed as the most reliable. It is worth mentioning that only recently QTAIM calculations with the use of basis functions (both of Slater or Gaussian type) of arbitrary angular momentum have been performed.17 However, this was done by means of a new Fortran code.17 In our studies Hartree-Fock (HF)14 and density functional theory (DFT)/B3LYP18-21 methods have been used. Results coming from this part of the calculations are, for better clarity, discussed separately for single, double, and triple covalent bonds in the next section. For formaldehyde we also performed QTAIM calculations for several values of CdO distances. In this case the B3LYP/ aug-cc-pVQZ optimal CdO distance of 1.1984 Å was changed by (0.01, (0.02, (0.03, (0.04, and (0.05 Å and kept constant during the geometry optimization while the rest of the geometrical parameters were reoptimized. QTAIM calculations with the CdO distance scan were then performed for four basis sets: 6-31++G(d,p), 6-311++G(d,p), aug-cc-pVDZ, and aug-ccpVTZ, together with the DFT/B3LYP approximation. Basis sets used in these calculations belong to two groups with respect to the description of the valence shell; namely, they are doubleor triple-ζ basis sets of either Pople or Dunning type. The aim of these calculations was to investigate the sensitivity of values of QTAIM parameters to a small change of the length of a covalent bond, since different values of bond length are what in fact one may find if the same molecule is settled in various chemical environments. The CdO bond seems to be a good object of interest for such analysis, since it is present in many molecular systems of large biological importance. It also may be considered as a good proton acceptor in the hydrogen bond phenomenon; thus as such it has an invaluable role in determining the structure of complex molecular systems. Besides, being double, it is between formally single and triple bonds. Results of these calculations are discussed in detail in a separate section. Calculations were performed with the use of the Gausian03 set of codes.22 A detailed analysis of the electron distribution function was made according to the concept of quantum theory of atoms in molecules1-6 using the AIMAll23 program. For acetylene, the use of some levels of approximation led to the improper characteristics of the CP of the CtC bond; namely, cage critical point (CCP) instead of BCP was obtained in these cases. Since in such a case AIMAll automatically interrupts calculation reporting unphysical electron density topology, the Xaim24 program was used to extract the electron density parameters of the triple CC bond. Properties of the electron density calculated at bond critical points of relevant covalent bonds were characterized. The electron density at BCP, FBCP, its Laplacian, ∇2FBCP, and the total electronic energy density, HBCP, were discussed. The list of these parameters, which are the most often used in studies of chemical interactions by means of QTAIM, was completed, with values of the delocalization index, δ(X,Y), which is the magnitude of the exchange of electrons in the basin of an atom X with those in the basin of an atom Y and is obtained by integration over X and Y atomic basins. Basis Set and Method Dependence. The aim of this paper is to analyze how values of the most important QTAIM parameters, being used in the description of a wide range of chemical interactions, behave if one changes the basis set and the method. Thus, we intend to investigate the reliability of

J. Phys. Chem. A, Vol. 114, No. 47, 2010 12499 QTAIM calculations in the theoretical description of covalent bonds. This “reliability” should be, however, understood as the stability of values of diverse QTAIM parameters with respect to the change of both the basis set and the method (HF or DFT/ B3LYP). For the purpose of our analysis we have chosen several simple molecules possessing formally single, double, or triple covalent bonds. Then parameters of the electron density function (F, ∇2F, H) were computed in BCP corresponding to a relevant covalent bond by means of diverse basis sets and either HF or DFT/B3LYP method. The delocalization index, δ(X,Y), was also calculated by integration over X and Y atomic basins. Case of Single Bonds. The selected group of molecules possessing single bonds consists of ethane (CsC), methanol (CsO), and fluoromethane (CsF). Simultaneously these bonds represent various extents of polarization due to the presence of a distinct element being the bonding partner of the carbon atom. Values of FBCP, its Laplacian (∇2FBCP), HBCP, and δ(C,Y) (Y ) C, O, F) obtained by means of various levels of approximation are shown in Figure 1. For more detailed plots showing separately dependences between QTAIM parameters and the level of approximation for relevant single bonds, see Figures S1-S4 in the Supporting Information. Analyzing dependences of the electron density in BCPs of CsC, CsO, and CsF single bonds in ethane, methanol, and fluoromethane, respectively, on the basis set used in QTAIM HF calculations, one can see that for ethane FB3LYP BCP (CC) < FBCP(CC), whereas this relation is opposite in the case of methanol and fluoromethane, thus for polarized covalent bonds. This finding indicates that HF method leads to the greater concentration of the electron density in the BCP of the C-C homonuclear single bond than DFT/B3LYP method does. As a consequence HF method describes this bond as being stronger than if QTAIM analysis were performed by means of DFT/B3LYP method. This conclusion has its confirmation in a more negative value of HF (CC) (see Figure 1, right-hand side, top). The opposite ∇2FBCP is true in the case of polarized C-O and C-F bonds in methanol and fluoromethane, respectively. Although, in general, a large dependence of FBCP on the basis set is noticeable, this dependence is similar for both HF and DFT/B3LYP methods, i.e., values of FBCP obtained by means of either HF or DFT/B3LYP follow the same trends. If one assumes that the largest aug-cc-pVTZ basis set gives the most reliable values of analyzed QTAIM parameters12 in BCPs of relevant covalent bonds, then one sees from Figure 1 that similar values of FBCP are also obtained if the largest 6-311++G(2df,2pd) (8) and 6-311++G(3df,3pd) (9) basis sets of Pople type are used. Also much smaller 6-31G(d,p) (3) basis set gives values of FBCP close to those from aug-cc-pVTZ. Just the opposite is found for the smallest 6-31G basis set (1) which in both HF and DFT/B3LYP gives considerably underestimated values of FBCP. Also poor results are obtained if 6-311++G(d,p) (7), ccpVDZ (10), and in particular aug-cc-pVDZ (11) basis sets are used (see Figure 1, left-hand side, top). The above findings are more or less common for all single bonds investigated in this study. Analyzing Figure 1 one can notice that, investigating for example FB3LYP BCP values obtained by means of aug-cc-pVTZ basis set, this value is the greatest for C-O. It may result from the partial sharing of electron density in the BCP of the C-O bond from the electron lone pairs of oxygen which is not so much concentrated as in the case of fluorine atom which is both smaller and more electronegative. Similarly as previously described, opposite dependences HF B3LYP and FBCP obtained for the homonuclear C-C between FBCP

12500

J. Phys. Chem. A, Vol. 114, No. 47, 2010

Jabłon´ski and Palusiak

Figure 1. Values of FBCP, ∇2FBCP, HBCP, and δ(C,Y) (Y ) C, O, F) for single bonds in ethane (CC), methanol (CO), and fluoromethane (CF) for HF and B3LYP and various basis sets. Basis sets are given numerical values: (1) 6-31G, (2) 6-31G(d), (3) 6-31G(d,p), (4) 6-31+G(d), (5) 6-31+G(d,p), (6) 6-31++G(d,p), (7) 6-311++G(d,p), (8) 6-311++G(2df,2pd), (9) 6-311++G(3df,3pd), (10) cc-pVDZ, (11) aug-cc-pVDZ, (12) cc-pVTZ, and (13) aug-cc-pVTZ. Points have been connected for a better visualization of results.

single bond and either C-O or C-F (see Figure 1, left-hand side, top), such opposite relations are also obtained in the case of ∇2FBCP (see Figure 1, right-hand side, top). However, for HF B3LYP (CC) < ∇2FBCP (CC) in the case Laplacian one obtains ∇2FBCP of the CsC bond in ethane, whereas for CsO in methanol and HF B3LYP > ∇2FBCP . The former CsF in fluoromethane one gets ∇2FBCP relation shows that HF method favors the covalent character of the CsC bond. Interestingly enough, in the case of CsF bond in fluoromethane, for most of the investigated basis sets, one obtains opposite signs of ∇2FBCP(CF) in HF and DFT/B3LYP (see Figure 1, right-hand side, top). The Laplacian of FBCP(CF) is negative in DFT/B3LYP, whereas it is positive in HF. Since, in general, ∇2FBCP is positive for close-shell interactions5,6 the positive value of ∇2FBCP obtained at the HF level indicates that HF method overestimates the ionic character of polarized covalent single bonds. This is also valid for C-O in methanol, since ∇2FBCP(CO) has considerably less negative value for HF calculations than for DFT/B3LYP. The dependence of ∇2FBCP on the basis set seems to be more complicated than it was found for FBCP. This is particularly veracious for polar bonds. This vague dependence is because of ∇2FBCP, since being the second derivative property, it is much more sensitive to small changes in the electron density distribution in the X-Y region. Nevertheless one can see from Figure 1 (right-hand side, top) that, at the HF level, 6-31G and ccpVDZ basis sets give values much different than if aug-cc-pVTZ basis set is used (this, however, is not valid for CsO, since HF (CO) very close to that of aug-cc6-31G gives value of ∇2FBCP pVTZ). Also the largest Pople-type basis sets 8 and 9 perform

not so well as in the case of FBCP if polarized single bonds are considered. At the DFT/B3LYP level this relationship is more diversified and depends in greater amount on the type of a single bond. For example 6-31G gives much less negative values of B3LYP (with respect to aug-cc-pVTZ) in the case of C-C ∇2FBCP and C-O bonds, but it performs well in the case of C-F bond. Also cc-pVDZ performs poorly for CsC and CsF, but it gives a value similar to that obtained by means of aug-cc-pVTZ for the C-O bond. Comparing values of ∇2FBCP obtained for CsC, CsO, and CsF single bonds in ethane, methanol, and fluoromethane, respectively, one sees that ∇2FBCP for CsC and CsO have similar values, whereas ∇2FBCP for CsF is either negative and close to zero (DFT/B3LYP) or positive (HF); thus the QTAIM analysis reveals the substantial ionic character of the CsF bond in fluoromethane. The total electronic energy density (H) belongs to the most important QTAIM parameters used in studies of chemical interactions. Again, as shown in Figure 1 (left-hand side, bottom), the usage of 6-31G, cc-pVDZ, and aug-cc-pVDZ basis sets leads to HBCP values which are clearly different than those obtained by means of the largest aug-cc-pVTZ. Also 6-311++ G(d,p) (7) basis set gives values similar to those obtained by means of DZ-quality Dunning-type basis sets (10 and 11), all being less negative than those from aug-cc-pVTZ. This finding is particularly valid for polarized CsO and CsF bonds. It seems that the enlargement of the set of polarization functions in the case of 6-311++G(d,p) (7), i.e. the transition to 8 or even better to 9, leads to HBCP values being closer to those from aug-ccpVTZ.

QTAIM Calculations for Covalent Bonds

J. Phys. Chem. A, Vol. 114, No. 47, 2010 12501

Figure 2. Values of FBCP, ∇2FBCP, HBCP, and δ(C,Y) (Y ) C, O) for double bonds in ethene (CC) and formaldehyde (CO) for HF and B3LYP and various basis sets. Basis sets are given numerical values: (1) 6-31G, (2) 6-31G(d), (3) 6-31G(d,p), (4) 6-31+G(d), (5) 6-31+G(d,p), (6) 6-31++G(d,p), (7) 6-311++G(d,p), (8) 6-311++G(2df,2pd), (9) 6-311++G(3df,3pd), (10) cc-pVDZ, (11) aug-cc-pVDZ, (12) cc-pVTZ, and (13) aug-cc-pVTZ. Points have been connected for a better visualization of results.

The delocalization index, δ(X,Y), provides a quantitative measure of the sharing of electrons between X and Y atomic basins. In other words, it provides a measure of the extent of delocalization of the electrons from one atomic basin into another.25 Thus it should be equal to 1 for an equally shared pair of electrons as it takes place in the case of homonuclear single bond. From the data shown in Figure 1 (right-hand side, bottom) it is seen that unity is almost approached for TZ-quality basis sets of Dunning type at the HF level. Similar value is also obtained if one computes δ(C,C) by means of 6-311++G(d,p) (7) basis set, whereas enlarging the set of polarization functions, i.e., transition to 6-311++G(2df,2pd) (8) or 6-311++G(3df,3pd) (9), leads to somewhat smaller value of δ(C,C) in ethane. The opposite trend is observed in the case of DFT/B3LYP calculations; the transition from 7 to 8 or 9 leads to δ(C,C) being closer to unity. It can be seen from Figure 1 that this time the worst δ(C,C) values are obtained if one employs 6-31G(d) (2) or 6-31+G(d) (4) basis sets. The addition of polarization functions for hydrogens significantly corrects these values. Considering results obtained for polarized C-O and C-F bonds, one may conclude that 6-31G basis set gives unsatisfactory values of δ(C,Y) (Y ) O, F). Similarly, cc-pVDZ and, toin a lesser extent, 6-311++G(d,p) lead to values which are considerably different than those obtained by means of augcc-pVTZ. In the case of 6-311++G(d,p) (7) this result can be corrected by the use of a double or triple set of polarization functions. Comparing δ(C,Y) values obtained by means of either 2 and 3 or 4 and 5, it seems that the presence of polarization functions on hydrogens has some importance. This is even more observable in the case of C-C homonuclear bond in ethane as

mentioned earlier. Finally it is clear from Figure 1 that all these findings are valid for both HF and DFT/B3LYP calculations due to similar trends of both with respect to the basis set. It is also seen that QTAIM calculations reproduce the decrease of the covalent character of investigated single bonds, since values of δ(C,Y) decrease in the CsC, CsO, and CsF order. Case of Double Bonds. We have investigated two types of double bonds. These are the homonuclear CdC bond in ethene and the heteronuclear CdO bond in formaldehyde. Results which show the influence of the basis set used in QTAIM calculations on values of FBCP, ∇2FBCP, HBCP, and δ(C,Y) (Y ) C, O) are presented in Figure 2 (see also more detailed Figures S5-S8 in the Supporting Information). For the CdO bond in formaldehyde one obtains almost the same values of FBCP(CO) irrespectively whether HF or DFT/ B3LYP QTAIM calculations are performed. For CdC in ethene DFT/B3LYP gives smaller values; however, in both cases the dependence of FBCP on the basis set is similar (see Figure 2, left-hand side, top). If we assume that aug-cc-pVTZ basis set gives the most reliable12 values of FBCP, then again poor results are obtained if QTAIM calculations are performed by means of 6-31G, cc-pVDZ and aug-cc-pVDZ basis sets. Comparing FBCP values obtained by means of cc-pVDZ and aug-cc-pVDZ basis sets, one sees that augmentation of the cc-pVDZ basis set leads to a small correction of FBCP values; i.e., it leads to FBCP values being closer to those obtained if the TZ-quality Dunningtype basis sets were used. This effect of augmentation of the cc-pVDZ basis set was not found in the case of single bonds (see Figure 1, left-hand side, top). Such a small correction of FBCP values due to the presence of diffuse functions is also

12502

J. Phys. Chem. A, Vol. 114, No. 47, 2010

observable if one compares 6-31G(d) (2) with 6-31+G(d) (4) or 6-31G(d,p) (3) with 6-31+G(d,p) (5). The best (i.e., the closest to those from aug-cc-pVTZ QTAIM calculations) FBCP values obtained by the use of a double-ζ-quality Pople-type basis set are for 6-31++G(d,p) (6) basis set. Interestingly enough, splitting the valence shell description to triple-ζ, i.e., the use of 6-311++G(d,p) (7), leads to worse FBCP values, which, however, are corrected if one uses the double (or triple) set of polarization functions (see basis sets 8 and 9 in Figure 2). Similarly as it was found for the CsC and CsO pair of single bonds in ethane and methanol, respectively, one obtains greater electron density values in the BCP of the CdO bond in formaldehyde than in the BCP of the CdC bond in ethene. This can be explained by the delocalization of the oxygen electron lone pairs in the CdO bond region. Similarly as it was found for FBCP, also in the case of ∇2FBCP, cc-pVDZ and aug-cc-pVDZ basis sets lead to poor results. Moreover, in the case of the CdC bond in ethene the smallest 6-31G basis set also tends to give a value of ∇2FBCP(CC) very close to those obtained by means of cc-pVDZ and aug-cc-pVDZ basis sets. However, this finding is not valid in the case of the heteronuclear CdO bond, where, particularly at the HF level, one obtains ∇2FBCP(CO) value practically the same as from TZquality Dunning-type basis sets. Similar finding applies to the C-O single bond in methanol (compare Figure 2 and Figure 1). Very close to those obtained by means of cc-pVTZ and augcc-pVTZ basis sets, values of ∇2FBCP(CO) are obtained if the 6-311++G(2df,2pd) basis set is used. In the case of HBCP parameter one obtains some interesting differences in the dependence on the basis set depending on whether CdC or CdO bond is investigated. Namely, while both 6-31G and DZ-quality Dunning-type basis sets give poor HBCP values, in the case of CdC, values obtained by means of DZquality Dunning-type basis sets are significantly less negative than those obtained by means of small and medium-sized Popletype basis sets (2-7), whereas in the case of CdO in formaldehyde relevant HBCP values are similar to each other (see Figure 2, left-hand side, bottom). Also in the case of the CdC bond in ethene, 6-311++G(2df,2pd) and to a somewhat lesser extent 6-311++G(3df,3pd) basis sets give HBCP values which are the closest to those obtained by means of 12 and 13 among all basis sets investigated in this study. This, however, is not the case for the CdO bond in formaldehyde. For these basis sets, in this case, one obtains significantly lower values of HBCP. If one compares HBCP values calculated for BCPs of CdC and CdO double bonds, then HBCP(CdO) < HBCP(CdC). The same relation is obtained for C-C and C-O single bonds (see Figure 1, left-hand side, bottom): i.e., HBCP(C-O) < HBCP(C-C). It means that for the same type of bond (single, double) the polarity of a bond leads to a lower value of HBCP. Conclusions applying to δ(C,Y) (Y ) C, O) values for investigated double bonds are similar to those already described for single bonds. Namely, for CdC in ethene, the worst results are obtained if one uses 6-31G(d) and 6-31+G(d) basis sets. However, it seems that the addition of the polarization functions for hydrogens corrects this results. For CdO in formaldehyde considerably poor result is obtained in the case in which the 6-31G basis set is used. Case of Triple Bonds. The group of systems possessing triple bonds is represented by three molecules. These are acetylene (HCtCH), nitrogen molecule (NtN), and hydrogen cyanide (HCtN). Dependences of FBCP, ∇2FBCP, HBCP, and δ(X,N) (X ) N, C) calculated for relevant triple bonds on the basis set used in QTAIM calculations are shown in Figure 3 (more

Jabłon´ski and Palusiak detailed plots can be found in Supporting Information as Figures S9-S12). For all, CtC, NtN, and CtN triple bonds investigated in this study, HF and DFT/B3LYP methods give very similar values of FBCP (see Figure 3, left-hand side, top). The worst values of FBCP (comparing to the aug-cc-pVTZ results) are obtained for the smallest 6-31G (1) basis set. In the case of the CtN and particularly CtC triple bond in hydrogen cyanide and acetylene, respectively, similarly poor FBCP value is obtained if the cc-pVDZ basis set (10) is used. Augmentation of the ccpVDZ basis set leads to the improvement of FBCP values. Again values of FBCP which are close to those if TZ-quality Dunningtype basis sets are used can be also obtained if the largest 8 or 9 Pople-type basis sets are taken for the QTAIM calculations. HF and DFT/B3LYP methods lead to similar values of ∇2FBCP if the triple bond is not polar (see Figure 3, right-hand side, top). In the case of the polar CtN bond in hydrogen cyanide, some more noticeable differences between HF and DFT/B3LYP values are seen. Particularly this is valid for the largest basis sets used in our investigations. If one compares ∇2FBCP values obtained by means of the largest TZ-quality Dunning-type basis sets with those obtained with the use of smaller basis sets, then it is clearly seen from Figure 3 that, again, poor results are obtained if either the cc-pVDZ or aug-cc-pVDZ basis set is used. Only 6-311++G(2df,2pd) and in a somewhat lesser extent 6-311++G(3df,3pd) basis sets give values which are close to those of 12 or 13. However, in the case of the polar CtN bond in hydrogen cyanide, also 6-311++G(d,p) and 6-31G basis sets give similar results. From Figure 3 it is also seen that the most negative value of ∇2FBCP is obtained for the NtN triple bond in the nitrogen molecule, thus indicating the largest covalent character among triple bonds investigated here. Interestingly enough, in the case of polarized CtN bond in the hydrogen cyanide one obtains much smaller values of ∇2FBCP which, depending on the basis set, can be close to zero or even positive. Such values of Laplacian of the electron density in the BCP of the CtN bond in the HCN molecule can be explained by the significant polarization of the CtN bond due to the presence of both highly electronegative nitrogen atom and π-type bonds. As a consequence this bond gains significantly ionic character and as such is described by close to zero or positive values of ∇2FBCP (CN) in the BCP of the CtN bond. For HBCP one obtains some differences in dependence on the basis set used in QTAIM calculations depending on whether or not the investigated triple bond is polar. While for CtC and NtN in acetylene and nitrogen molecule, respectively, 1, 10, and 11 basis sets give poor results (as compared to aug-ccpVTZ), in the case of CtN in HCN, DZ-quality basis sets of Dunning-type perform better than the smallest Pople-type basis sets (1-6). It is also seen that 6-311++G(2dp,2pd) basis set (8) leads to values of HBCP which are close to those of aug-ccpVTZ. HF and DFT/B3LYP methods also show (see Figure 3, righthand side, bottom) similar dependence on the basis set in the case of δ(X,N) values. If one considers the NtN bond in the nitrogen molecule then, again, results obtained by means of ccpVDZ basis set are somewhat distinct (see also Figure S12) from those of aug-cc-pVTZ, whereas the aug-cc-pVDZ basis set performs pretty well. This is not, however, the case for the CtN polar bond in the HCN molecule. Now both cc-pVDZ and aug-cc-pVDZ basis sets give poor results, as compared to aug-cc-pVTZ, and only 7-9 basis sets give δ(C,N) values being close to those from aug-cc-pVTZ (or cc-pVTZ).

QTAIM Calculations for Covalent Bonds

J. Phys. Chem. A, Vol. 114, No. 47, 2010 12503

Figure 3. Values of FBCP, ∇2FBCP, HBCP, and δ(X,N) (X ) N, C) for triple bonds in acetylene (CC), nitrogen molecule (NN), and hydrogen cyanide (CN) for HF and B3LYP and various basis sets. Basis sets are given numerical values: (1) 6-31G, (2) 6-31G(d), (3) 6-31G(d,p), (4) 6-31+G(d), (5) 6-31+G(d,p), (6) 6-31++G(d,p), (7) 6-311++G(d,p), (8) 6-311++G(2df,2pd), (9) 6-311++G(3df,3pd), (10) cc-pVDZ, (11) aug-cc-pVDZ, (12) cc-pVTZ, and (13) aug-cc-pVTZ. Points have been connected for a better visualization of results.

Influence of the Bond Length. The most noticeable dependence between FBCP, ∇2FBCP, HBCP, δ(X,Y), and the basis set used in QTAIM calculations has been found in the case of both multiple and polar bonds. Thus, the CdO bond in formaldehyde seems to be a proper structural unit in further studies of the influence of the small change in the bond length to the value of the relevant QTAIM parameter. Such a study is important, since, in fact, a given bond does not have a unique length, but its length depends on the chemical environment due to the different perturbation in the electron density distribution in the bond region. Thus, the CdO double bond in formaldehyde is a good choice also from this point of view, since as an effective proton acceptor in the hydrogen bond, it has various lengths in different chemical situations. Moreover it is an important building block of many molecular systems (consider, for example, the carboxylic or peptide group). To study the influence of the small perturbation in the CdO bond distance we have performed scan calculations for various values of CdO bond length as described in Methodology. To mention briefly, the B3LYP/aug-cc-pVQZ optimal distance of 1.1984 Å was changed by (0.01, (0.02, (0.03, (0.04, and (0.05 Å and kept constant during the geometry optimization procedure. The rest of the geometrical parameters were reoptimized on the B3LYP/aug-cc-pVQZ level of theory. Then, QTAIM calculations were performed for four basis sets: 6-31++G(d,p), 6-311++G(d,p), aug-cc-pVDZ, and aug-ccpVTZ by means of the DFT/B3LYP method. Results are shown in Figure 4. As it is seen from Figure 4, the compression of the CdO bond in formaldehyde leads to the increase of the electron

density, its Laplacian, and the absolute value of the total electronic energy density, while the delocalization index decreases. This picture is obtained for all basis sets. These changes, being a consequence of the CdO bond compression, indicate the larger amount of the electron density concentration in the BCP of the CdO bond, but simultaneously the larger separation of the charge between C and O attractors what, as a consequence, explains the increase of ∇2FBCP(CdO) values. Also it explains the decrease of the δ(C,O) value (see Figure 4, righthand side, bottom). Although the decrease of the delocalization index upon the CdO bond compression may seem to be somewhat surprising, the increased separation of charge in the CdO bond may be explained by the stronger overlapping of orbitals of C and O atoms, mainly of p orbitals forming the π-bond, thus leading to the facilitated charge flow toward the more electronegative oxygen atom. In other words, the increasing ionicity of the CdO bond decreases the delocalization index due to unequal sharing between the two basins. Results shown in Figure 4 also indicate that point distributions obtained by means of 6-31++G(d,p) and aug-cc-pVDZ are close to each other, while, at least in the case of ∇2FBCP(CdO) and HBCP(CdO), point distributions obtained by means of augcc-pVTZ are closer to those obtained by means of the 6-311++G(d,p) basis set. This is understood since in both cases relevant basis sets belong to the same group with respect to the description of the valence shell; namely, in the former case they belong to the group of split-valence double-ζ basis sets, while in the latter they belong to the group of split-valence triple-ζ basis sets.

12504

J. Phys. Chem. A, Vol. 114, No. 47, 2010

Jabłon´ski and Palusiak

Figure 4. Values of FBCP, ∇2FBCP, HBCP, and δ(C,O) obtained for various distances of CdO double bond in formaldehyde calculated by means of B3LYP method and several basis sets. Points have been connected for a better visualization of results.

It should be mentioned that very similar point distributions indicating the dependence between values of relevant QTAIM parameters and the CdO bond length were obtained when geometry optimizations of formaldehyde were performed by means of the same basis sets for which also QTAIM analysis was done (see Figure S17 in Supporting Information). Although this way of proceeding could be considered more justified due to the fact that exactly the same method and the basis set are used for both the geometry optimization and QTAIM calculations, it is not in agreement with the assumption of the present calculations where we studied the influence of a basis set to values of QTAIM parameters if the structure of formaldehyde is fixed and thus exactly the same for all basis sets. It is noteworthy that even for such small changes of the CdO bond length ∇2FBCP(CdO) in the BCP of this bond may accept values from a rather wide range. For example, taking into account values of ∇2FBCP(CdO) obtained by means of B3LYP/ aug-cc-pVTZ calculations, one may encounter negative as well as positive values of ∇2FBCP(CdO) and the spread of ∇2FBCP(CdO) values is about 0.7 au. This spread is even larger in the case of the aug-cc-pVDZ basis set and amounts to 1.4 au for the investigated range of changes in the CdO bond distances. Thus, taking into account that a given structural unit (e.g., CdO bond) can be in diverse chemical situations, it seems to be deprived of grounds to discuss about the correct value of Laplacian of the electron density in the BCP of the CdO bond in many molecular systems.26 Conclusions The influence of the basis set and the level of theory to values of FBCP, ∇2FBCP, HBCP, and δ(X,Y) was studied by computing

relevant QTAIM parameters for several covalent bonds for a large set of small- and medium-sized basis sets of Pople and Dunning type and by means of HF and DFT/B3LYP methods. The following most important conclusions can be formulated on the basis of our analysis: (a) In general, HF and DFT/B3LYP methods give very similar distributions of values of relevant QTAIM parameters with respect to the basis set. (b) In many cases, for the same basis set, HF and DFT/ B3LYP methods lead to very similar values of FBCP, ∇2FBCP, and HBCP. This finding in particular applies to triple bonds. (c) The discrepancy of values of relevant QTAIM parameters with respect to the basis set is the most noticeable for double and triple bonds, while for single bonds this discrepancy is much smaller. (d) The discrepancy of values of relevant QTAIM parameters with respect to the basis set is more noticeable for polar bonds than for nonpolar bonds. (e) The largest influence of the basis set is obtained for the Laplacian of the electron density (∇2FBCP). (f) In general, the worst results are obtained by means of the 6-31G basis set and to a somewhat lesser extent cc-pVDZ and aug-cc-pVDZ basis sets. (g) In general, values of relevant QTAIM parameters being the closest to those obtained by means of either cc-pVTZ or aug-cc-pVTZ basis set were also obtained if 6-311++G(2df,2pd) and 6-311++G(3df,3pd) basis sets were used. Moreover 6-311++G(2df,2pd) seems to perform better than 6-311++G(3df,3pd). This finding in particular applies to double and triple CC bonds and for ∇2FBCP and HBCP parameters.

QTAIM Calculations for Covalent Bonds (h) Values of QTAIM parameters reflect the increasing strengths of bonds as follows: CsC < CdC < CtC and CsO < CdO. The findings listed above show that QTAIM analysis of chemical interactions should not be performed by means of 6-31G and DZquality Dunning-type basis sets, i.e., neither cc-pVDZ nor aug-ccpVDZ. On the contrary, we recommend aug-cc-pVTZ basis set or either 6-311++G(2df,2pd) or 6-311++G(3df,3pd); however, the smaller 6-311++G(d,p) basis set also may be recommended if the use of larger basis sets is too time-consuming. The use of ccpVDZ and aug-cc-pVDZ basis sets has also not been recently recommended in studies of weak hydrogen bonds by means of QTAIM calculations. Present results also show that QTAIM calculations seem to be much more method/basis set dependent in the case of covalent bonds than it was found in the case of weak hydrogen bonds. Finally, it was shown that the compression of the CdO bond in formaldehyde leads to the increase of the electron density, its Laplacian, and the absolute value of the total electronic energy density, while the delocalization index decreases. These changes of relevant QTAIM parameters can be explained by the increasing ionic character of the CdO bond during compression. The same picture was obtained for 6-31++G(d,p), 6-311++G(d,p), aug-cc-pVDZ, and aug-cc-pVTZ basis sets and by means of DFT/B3LYP method. Acknowledgment. Calculations have been carried out on the multiprocessor cluster at the Information and Communication Technology Center of Nicolaus Copernicus University, Torun´ (http://www.ucu.umk.pl). M.P. acknowledges University Grant No. 505/721/R. Supporting Information Available: Detailed figures presenting dependences between values of FBCP, ∇2FBCP, HBCP, and δ(X,Y) calculated by means of relevant basis sets and by means of either HF or DFT/B3LYP method for single (Figures S1-S4), double (Figures S5-S8), and triple (Figures S9-S12) bonds, for single (ethane), double (ethene), and triple (acetylene) carbon-carbon bonds (Figure S13), and for single (methanol) and double (formaldehyde) carbon-oxygen bonds (Figure S14; both Figures S13 and S14 present the influence of the bond order to the dependence between values of relevant QTAIM parameters and both the basis set and the method used in QTAIM calculations), showing either C-C (ethane) and C-O (methanol) single bonds or for CdC (ethene) and CdO (formaldehyde) double bonds (Figures S15 and S16, respectively; these figures facilitate analysis of the influence of the

J. Phys. Chem. A, Vol. 114, No. 47, 2010 12505 bond polarity to the distribution of values of relevant QTAIM parameters with respect to the level of approximation used in QTAIM calculations), and depicting the result of scan calculations for formaldehyde if both the geometry optimization and the QTAIM calculations were performed on the same level of theory (Figure S17; see the text for explanation of this approach), and tables listing geometries of ethane, fluoromethane, methanol, ethene, formaldehyde, acetylene, nitrogen molecule, and hydrogen cyanide obtained at the B3LYP/aug-cc-pVQZ level of theory and geometries of formaldehyde obtained by means of DFT/B3LYP method and 6-31++G(d,p), 6-311++G(d,p), augcc-pVDZ, and aug-cc-pVTZ basis sets. This information is available free of charge via the Internet at http://pubs.acs.org. References and Notes (1) Bader, R. F. W. Acc. Chem. Res. 1975, 8, 34. (2) Bader, R. F. W.; McDougall, P. J.; Lau, C. D. H. J. Am. Chem. Soc. 1984, 106, 1594. (3) Bader, R. F. W.; Esse´n, H. J. Chem. Phys. 1984, 80, 1943. (4) Bader, R. F. W. Acc. Chem. Res. 1985, 18, 9. (5) Bader, R. F. W. Atoms in Molecules: A Quantum Theory; Oxford University Press: New York, 1990. (6) Popelier, P. L. A. Atoms in Molecules. An Introduction; Longman: Singapore, 2000. (7) Knop, O.; Boyd, R. J.; Choi, S. C. J. Am. Chem. Soc. 1988, 110, 7299. (8) Knop, O.; Rankin, K. N.; Boyd, R. J. J. Phys. Chem. A 2001, 105, 6552. (9) Knop, O.; Rankin, K. N.; Boyd, R. J. J. Phys. Chem. A 2003, 107, 272. (10) Castillo, N.; Robertson, K. N.; Choi, S. C.; Boyd, R. J.; Knop, O. J. Comput. Chem. 2008, 29, 367. (11) Matta, C. F. J. Comput. Chem. 2010, 31, 1297. (12) Jabłonski, M.; Palusiak, M. J. Phys. Chem. A 2010, 114, 2240. (13) Mata, I.; Alkorta, I.; Molins, E.; Espinosa, E. Chem.sEur. J. 2010, 16, 2442. (14) Szabo, A.; Ostlund, N. S. Modern Quantum Chemistry: Introduction to AdVanced Electronic Structure Theory; Macmillan: New York, 1982. (15) Kendall, R. A.; Dunning, T. H., Jr.; Harrison, R. J. J. Chem. Phys. 1992, 96, 9747. (16) Davidson, E. R. Chem. Phys. Lett. 1996, 260, 514. (17) Volkov, A.; Koritsanszky, T.; Chodkiewicz, M.; King, H. F. J. Comput. Chem. 2009, 30, 1379. (18) Hohenberg, P.; Kohn, W. Phys. ReV. 1964, 136, B864. (19) Kohn, W.; Sham, L. J. Phys. ReV. 1965, 140, A1133. (20) Becke, A. D. J. Chem. Phys. 1993, 98, 5648. (21) Lee, C.; Yang, W.; Parr, R. G. Phys. ReV. B 1988, 37, 785. (22) Frisch, M. J. et. al.; Gaussian 03; Gaussian: Wallingford, CT, 2004. (23) Keith, T. A. AIMAll (Version 10.06.21), aim.tkgristmill.com, 2010. (24) Ortiz, C. O.; Bo, C. Xaim, Universitat Rovira i Virgili, Tarragona. Spain, http://www.quimica.urv.es/XAIM/, accessed November 2009. (25) Fradera, X.; Austen, M. A.; Bader, R. F. W. J. Phys. Chem. A 1999, 103, 304. (26) Birkedal, H.; Madsen, D.; Mathiesen, R. H.; Knudsen, K.; Weber, H.-P.; Pattison, P.; Schwarzenbach, D. Acta Crystallogr. 2004, A60, 371.

JP106740E