ARTICLE pubs.acs.org/JPCA
Subshell Fitting of Relativistic Atomic Core Electron Densities for Use in QTAIM Analyses of ECP-Based Wave Functions Todd A. Keith* Semichem, Inc., 12456 West 62nd Terrace, Suite D Shawnee, Kansas 66216, United States
Michael J. Frisch Gaussian, Inc., 340 Quinnipiac Street, Building 40 Wallingford, Connecticut 06492, United States
bS Supporting Information ABSTRACT: Scalar-relativistic, all-electron density functional theory (DFT) calculations were done for free, neutral atoms of all elements of the periodic table using the universal Gaussian basis set. Each core, closed-subshell contribution to a total atomic electron density distribution was separately fitted to a spherical electron density function: a linear combination of s-type Gaussian functions. The resulting core subshell electron densities are useful for systematically and compactly approximating total core electron densities of atoms in molecules, for any atomic core defined in terms of closed subshells. When used to augment the electron density from a wave function based on a calculation using effective core potentials (ECPs) in the Hamiltonian, the atomic core electron densities are sufficient to restore the otherwise-absent electron density maxima at the nuclear positions and eliminate spurious critical points in the neighborhood of the atom, thus enabling quantum theory of atoms in molecules (QTAIM) analyses to be done in the neighborhoods of atoms for which ECPs were used. Comparison of results from QTAIM analyses with all-electron, relativistic and nonrelativistic molecular wave functions validates the use of the atomic core electron densities for augmenting electron densities from ECP-based wave functions. For an atom in a molecule for which a small-core or mediumcore ECPs is used, simply representing the core using a simplistic, tightly localized electron density function is actually sufficient to obtain a correct electron density topology and perform QTAIM analyses to obtain at least semiquantitatively meaningful results, but this is often not true when a large-core ECP is used. Comparison of QTAIM results from augmenting ECP-based molecular wave functions with the realistic atomic core electron densities presented here versus augmenting with the limiting case of tight core densities may be useful for diagnosing the reliability of large-core ECP models in particular cases. For molecules containing atoms of any elements of the periodic table, the production of extended wave function files that include the appropriate atomic core densities for ECP-based calculations, and the use of these wave functions for QTAIM analyses, has been automated.
’ INTRODUCTION Effective quantum chemical modeling of molecules or extended systems containing heavy atoms sometimes requires accounting for relativistic effects from the fast moving electrons of the heavy atom cores, the importance of these relativistic effects increasing rapidly with atomic number.13 Advances in relativistic quantum chemistry calculation methodologies and implementations, together with advances in hardware, have led to an increasing number of quantum chemical studies of heavy atom systems for which all electrons are treated explicitly and relativistically. Nevertheless, compared to nonrelativistic quantum chemistry calculations, relativistic calculations can still be expensive and limited in the size of systems one can study and in the properties one can calculate, due, for example, to difficulties in implementing and calculating analytic derivatives. r 2011 American Chemical Society
Effective core potentials (ECPs) derived from relativistic free atom calculations are commonly used in formally nonrelativistic quantum chemistry calculations of molecules and extended systems containing heavy atoms as a relatively straightforward, reliable and cost-effective way to implicitly account for some or all of the core electrons of the heavy atoms.3,4 Computational advantages of the ECP approach are fewer basis functions and fewer molecular orbitals, simpler energy functionals and more availability of analytic derivatives
Special Issue: Richard F. W. Bader Festschrift Received: April 29, 2011 Revised: June 30, 2011 Published: July 22, 2011 12879
dx.doi.org/10.1021/jp2040086 | J. Phys. Chem. A 2011, 115, 12879–12894
The Journal of Physical Chemistry A for optimizing and characterizing geometries and calculating molecular properties. The avoidance of explicitly treating some or all core electrons in calculations employing ECPs does have some drawbacks, however, beyond the obvious one of possibly reduced accuracy compared to an all-electron treatment. The lack of explicit representation of core electrons hinders the calculation and analysis of properties that depend on local or regional properties, such as the electron density distribution in the atomic core region. Important among these are properties determined using quantum theory of atoms in molecules (QTAIM) of Bader and co-workers.516 QTAIM analyses are based on the topological properties of electron density distributions, the dominant feature of which is strong maxima at the nuclear positions. This dominant topological feature leads to a unique and exhaustive partitioning of molecules and extended systems into basins of the electron density distribution’s gradient vector field, basins whose form, properties, and topological connectivity are consistent with the fundamental chemical ideas of atoms in molecules, functional groups and bonding. The absence of core electrons for an atom results in the absence of a corresponding nuclear maximum of the electron density distribution and therefore a qualitative change in the topology of the electron density distribution in the neighborhood of the atom, which makes meaningful QTAIM analyses difficult or impossible. Several studies1821 have shown that simply augmenting the electron density distribution obtained from a wave function of a Hamiltonian including ECPs (or even from semiempirical calculations17) with core electron densities obtained from allelectron, free atom calculations results in a correct topology of the molecular electron density distribution and thus enables one to perform QTAIM analyses. Comparison with results from allelectron calculations for a few small molecules has indicated that this procedure is also at least semiquantitatively reasonable.1821 Because ECPs are typically designed to model a core that includes relativistic effects,3,4 it is appropriate to use atomic core electron densities that are themselves obtained from relativistic all-electron calculations. Systematic calculations and tabulations of relativistic free atom core electron densities for all atoms of the periodic table for use in augmenting electron densities from ECP-based molecular calculations were presented by Cioslowski et al.18,19 Their approach, while quite useful, was incomplete in that only cores corresponding to noble gas electronic configurations were considered. This makes their data set unusable with a significant number of ECPs. For example, it is quite common for large-core ECPs of main group elements in the fourth row to include the 3d subshell electrons in the core. The electron density for such a core electron configuration, namely (1s2,2s2,2p6,3s2,3p6,3d10), is not available from the data set of Cioslowski et al. because it does not correspond to a noble gas configuration. Simple model core densities primarily for use in stabilizing density functional theory (DFT) calculations using certain ECPs have been given by Louie et al.22 and by Delley.23 It is the objective of the work described here to enable and, using QTAIM analyses, demonstrate results from a reliable and automated procedure for augmenting electron densities obtained from a broad range of ECP-based wave functions with appropriate relativistic atomic core electron densities.
ARTICLE
’ CALCULATION OF ELECTRON DENSITIES OF FREE ATOMS Large basis set, scalar-relativistic, all-electron KohnSham spin density functional theory calculations were done for free, neutral, ground state atoms for several spin states of all elements of the periodic table using a development version of the Gaussian quantum chemistry package.24 The density functional used was B3LYP.2528 The basis set used was the universal Gaussian basis set (UGBS).29 Scalar-relativistic effects were accounted for by using the DouglasKrollHess second-order scalar-relativistic Hamiltonian (DKH2).3033 The spin orbitals were unrestricted. The single-determinant, implicitly correlated B3LYP DFT method was used instead of an ab initio multideterminant correlated method such as CI because of the desire to obtain atomic orbital contributions to the electron density that could all be readily assigned to atomic electronic subshells. In most cases, the lowest energy spin state found for each element was selected for further use in isolating and fitting closed-subshell atomic electron densities. The few exceptions were low-lying spin states with spherically symmetric total electron density distributions that were considered more likely to resemble those in a molecular environment. For example, for chromium the heptet state corresponding to the valence electronic configuration (4s1,3d5) was used rather than the ground quintet state corresponding to (4s2,3d4). The final spin multiplicities used for each element are provided as Supporting Information in Table S1. The total electron density F(r;ZA) for an atom with atomic number ZA is given in terms of subshell density contributions Fn,l(r;ZA) as Fðr;ZA Þ ¼
n1
∑ ∑ Fn, l ðr;ZA Þ n¼1 l¼0
ð1Þ
where Fn,l(r;ZA) is the electron density contribution from atomic orbitals corresponding to the lth subshell of the nth principal shell. Because the wave functions are single-determinants, the subshell atomic orbitals and therefore the corresponding subshell electron densities are readily determined.
’ FITTING OF SUBSHELL CONTRIBUTIONS TO ELECTRON DENSITIES OF FREE ATOMS For each closed atomic subshell, the corresponding set of atomic orbitals was identified and the total electron density Fn,l(r;ZA) from these subshell atomic orbitals was fitted to a compact spherical density function Fn,l(r;ZA), defined as a linear combination of s-type primitive Gaussian functions ϕ centered on the nuclear position RA of the atom: F_ n, l ðr;ZA Þ ¼
∑i Ci ðn;l;ZA Þ ϕi ðr R A Þ
ð2Þ
For example, for a bromine atom, six core closed-subshell densities were fitted to the corresponding subshell densities of the bromine ground state doublet: F 1,0 (r;Br), F 2,0 (r; Br), F 2,1 (r; Br), F 3,0 (r; Br), F 3,1 (r; Br), and F 3,2 (r; Br) corresponding to the 1s2, 2s2, 2p6, 3s2, 3p6, and 3d10 core subshells. The set of s-type Gaussian functions ϕ used in the expansion of Fn,l(r;ZA) was derived from the UGBS for the element with 12880
dx.doi.org/10.1021/jp2040086 |J. Phys. Chem. A 2011, 115, 12879–12894
The Journal of Physical Chemistry A
ARTICLE
Table 1. Explanation of QTAIM Atomic and Bonding Property Symbolsa,b,c q(A)
net charge of atom A
N(A)
electron population of atom A
Vol(A)
volume enclosed by the interatomic surfaces of atom A and the 0.001 electron density isosurface
μintra(A)
intra-atomic dipole moment of atom A, i.e., the dipole moment of the atomic charge distribution, with the nucleus as origin
Q(A) LI(A)
traceless quadrupole moment tensor of atom A, e.g., Qxx(A) = 3Æ(x XA)2æA Æ(r RA)2æA localization index of atom A, the average number of electrons localized in atom A; by definition, includes all ECP-modeled core electrons of atom A
DI(A,B)
delocalization index (bond order) between atoms A and B, i.e., the average number of electrons delocalized (shared) between atoms A and B; by definition, does not include contribution from core electrons modeled with ECPs
K(A)
electronic kinetic energy of atom A
Nleak(A)
number of ECP-modeled core electrons of atom A minus the integral of the corresponding augmenting core electron density over basin of atom A
Fb(A|B)
electron density at bond critical point (BCP) between atoms A and B
Fr
electron density at a ring critical point (RCP)
Fc r2Fb(A|B)
electron density at a cage critical point (CCP) Laplacian of electron density at a BCP between atoms A and B
r2Fr
Laplacian of electron density at a RCP
r2Fc
Laplacian of electron density at a CCP
εb(A|B)
ellipticity of the BCP between atoms A and B. εb(A|B) = λ1/λ2 1 where λ1 and λ2 are the two negative principal curvatures (eigenvalues) of the Hessian of the electron density, with λ1 being the most negative one.
rb(A|B)
Euclidean distance from the nucleus of atom A to the BCP between atoms A and B
Rb(A|B)
bond path length between bonded atoms A and B
Re(A,B) Kb(A|B)
Euclidean distance between nuclei of atoms A and B Hamiltonian electronic kinetic energy density at BCP between atoms A and B
Kr
Hamiltonian electronic kinetic energy density at a RCP
Kc
Hamiltonian electronic kinetic energy density at a CCP
Gb(A|B)
Lagrangian electronic kinetic energy density at BCP between atoms A and B
Gr
Lagrangian electronic kinetic energy density at a RCP
Gc
Lagrangian electronic kinetic energy density at a CCP
V b(A|B)
virial field at BCP between atoms A and B
Vr Vc
virial field at a RCP virial field at a CCP
Fb(A|B)
Ehrenfest force density at BCP between atoms A and B
a
See refs 516. b In this paper, all QTAIM and other quantities are expressed in atomic units. c Properties such as K(A), Gb(A|B), and Fb(A|B) that are not determined by the electron density distribution do not include contributions from core electrons modeled with ECPs or their corresponding augmentative core densities.
atomic number ZA according to the prescription given by Yang, Rendell, and Frisch.34 The fitting procedure used was similar to the method used to generate density fitting basis sets for speeding up calculations of the electron repulsion J matrix in DFT calculations,3537 but with the essential additional requirement that the electron density distribution be non-negative: F_ n, l ðr;ZA Þ g 0
" r ∈ R3
Z
Z dr1
Z
dξ½n, l;ZA =dCi ¼ 2
dr2 ðF_ n, l ðr1 ;ZA Þ Fn, l ðr1 ;ZA ÞðF _ n, l ðr2 ;ZA Þ
Z Fn, l ðr2 ;ZA ÞÞjr1 r2 j1 þ λ½ dr F _ n, l ðr;ZA Þ N n, l
ð4Þ where Fn,l(r;ZA) is constrained to integrate to the number of electrons, Nn,l, of the n,l subshell.
Z
dr1
dr2 ϕi ðr RÞðF_ n, l ðr2 ;ZA Þ
Fn, l ðr2 ;ZA ÞÞjr1 r2 j1 þ λ
Z dr ϕi ðr R A Þ ¼ 0
" Ci
ð5Þ
ð3Þ
In the common “Coulomb metric” variant of the density fitting procedure used in this work,3537 the following Coulombic fitting error functional ξ[n,l;ZA] is minimized with the respect to the fitting coefficients C: ξ½n;l;ZA ¼
Minimization of ξ[n,l;ZA] with respect to each of the fitting coefficients gives the following necessary conditions:
Equation 5 defines a set of linear equations that can be solved for the fitting coefficients C using standard methods. Because the fitting coefficients C are not constrained to be non-negative, it is possible for the fitted subshell density obtained from solving eq 5 to be negative at some point in space. Several different schemes were therefore employed to ensure positive definiteness of Fn,l(r;ZA), eq 3. The final fitted subshell densities were taken from the scheme that gave the best non-negative fit for that subshell. All of these schemes involved evaluating F n,l(r;Z A) on a dense radial mesh from the nucleus out to a distance where the density was negligibly small. For example, in one approach, if F n,l(r;Z A) was negative at any point, then the primitive function making the largest negative contribution to the 12881
dx.doi.org/10.1021/jp2040086 |J. Phys. Chem. A 2011, 115, 12879–12894
The Journal of Physical Chemistry A
ARTICLE
Figure 1. Molecular graph, atomic charges and electron density at bond critical points of GeAdamantane using six different methods. Contours of the electron density in a σd symmetry plane are also shown. Small green spheres are bond critical points (BCPs). Small red spheres are ring critical points (RCPs). Small blue spheres are cage critical points (CCPs). Small pink spheres are non-nuclear attractor critical points (NNACPs). Paths originating at BCPs and terminating at nuclei or NNACPs are bond paths. (a) All-electron nonrelativistic wave function. (b) All-electron relativistic wave function. (c) Wave function from a small-core, 10-electron ECP calculation. (d) Wave function from a large-core, 28-electron ECP calculation. (e) Wave function from a small-core, 10-electron ECP calculation augmented with Ge atomic small-core electron densities. (f) Wave function from a large-core, 28-electron ECP calculation augmented with Ge atomic large-core electron densities.
density at the most negative point was discarded and the basic fitting procedure defined by eqs 2, 4, and 5 was repeated. For about half of all the subshells of all the elements, this
produced the best fit. Other strategies involved selecting the primitive to delete in various other ways, which are detailed in the Supporting Information Document S2. Each 12882
dx.doi.org/10.1021/jp2040086 |J. Phys. Chem. A 2011, 115, 12879–12894
The Journal of Physical Chemistry A
ARTICLE
Table 2. QTAIM Atomic and Bonding Properties of TiO2 Using Different Methodsa property
A
B
C
D
E
F
G
H
q(Ti)
1.854
1.865
1.847
1.542
1.847
2.025
1.847
2.025
q(O)
0.927
0.932
0.923
0.771
0.923
1.012
0.923
1.012
Vol(Mol)
405.3
405.8
402.4
399.2
402.4
394.9
402.4
399.2
Vol(Ti)
97.1
97.1
98.2
105.6
98.2
72.0
98.1
72.1
Vol(O)
154.1
154.3
152.1
146.8
152.1
161.5
152.1
161.5
μ(Mol)z
2.563
2.571
2.543
2.535
2.543
2.535
2.543
2.534
μintra(Ti)z
0.332
0.333
0.411
0.010
0.411
0.306
0.411
0.305
μintra(O)y μintra(O)z
0.010 0.159
0.013 0.164
0.004 0.124
0.068 0.074
0.004 0.124
0.285 0.335
0.004 0.124
0.285 0.335
Qxx(Ti)
0.070
0.112
0.136
0.527
0.136
0.602
0.129
0.601
Qyy(Ti)
0.523
0.515
0.531
0.258
0.532
0.312
0.537
0.313
Qzz(Ti)
0.593
0.627
0.668
0.782
0.668
0.914
0.666
0.914 1.185
Qxx(O)
1.334
1.341
1.480
1.559
1.480
1.185
1.480
Qyy(O)
1.726
1.737
1.790
1.910
1.790
1.601
1.790
1.601
Qzz(O)
0.393
0.396
0.310
0.350
0.310
0.416
0.310
0.416
LI(Ti) LI(O)
18.325 7.931
18.320 7.939
18.320 7.912
18.694 7.677
18.320 7.912
18.433 8.093
18.320 7.912
18.433 8.093
DI(Ti|O)
1.822
1.815
1.833
1.889
1.833
1.542
1.833
1.542
DI(O|O)
0.170
0.171
0.189
0.172
0.189
0.295
0.189
0.295
Nleak(Ti)
na
na
0.0000
0.1257
0.0000
0.0000
0.0000
0.0000
Fb(Ti|O)
0.2500
0.2494
0.2380
0.2738
0.2380
0.1806
0.2380
0.1806
r2Fb(Ti|O)
0.9186
0.9297
1.1366
0.7551
1.1350
0.1420
1.1350
0.1420
εb(Ti|O)
0.0283
0.0315
0.0203
0.0035
0.0203
0.0159
0.0203
0.0159
rb(Ti|O) rb(O|Ti)
1.6034 1.4865
1.6025 1.4875
1.5990 1.4910
1.6627 1.4273
1.5989 1.4911
1.4003 1.6918
1.5989 1.4911
1.4003 1.6918
Kb(Ti|O)
0.1510
0.1501
0.1006
0.1430
0.1006
0.1709
0.1006
0.1709
Gb(Ti|O)
0.3806
0.3825
0.3843
0.1744
0.3843
0.1354
0.3843
0.1354
V b (Ti|O)
0.5316
0.5326
0.4849
0.3174
0.4849
0.3062
0.4849
0.3062
|FbE(Ti|O)|
0.1120
0.1158
0.0498
0.2276
0.0500
0.5456
0.0500
0.5456
a
The geometry for all methods is the equilibrium geometry at B3LYP/cc-pVDZ, all-electron, nonrelativistic. A. B3LYP/cc-pVDZ, all-electron, nonrelativistic. B. B3LYP/cc-pVDZ, all-electron, relativistic (DKH2). C. B3LYP/cc-pVDZ[O],LANL2DZ[Ti], 10-electron ECP for Ti, 10-electron atomic core density on Ti. D. B3LYP/cc-pVDZ[O],LANL1DZ[Ti], 18-electron ECP for Ti, 18-electron atomic core density on Ti. E. B3LYP/ccpVDZ[O],LANL2DZ Ti], 10-electron ECP for Ti, 10-electron tight core density on Ti. F. B3LYP/cc-pVDZ[O],LANL1DZ Ti], 18-electron ECP for Ti, 18-electron tight core density on Ti. G. B3LYP/cc-pVDZ[O],LANL2DZ Ti], 10-electron ECP for Ti, no core density on Ti. H. B3LYP/cc-pVDZ[O], LANL1DZ Ti], 18-electron ECP for Ti, no core density on Ti.
of the other eight strategies was found to be optimal for at least a few of the subshells of some elements
’ CONSTRUCTION OF CORE DENSITIES OF ATOMS IN MOLECULES Once all of the subshell density functions F n,l (r;Z A ) are obtained for an atom, it is straightforward to combine them to obtain the total core electron density F core(r;Z A ,N core,A ) for any N core,A electron core for the atom that is defined in terms of closed subshells. For example, a 28 electron largecore density F c(r;Se,28) for selenium atoms, corresponding to electron configuration (1s2,2s2,2p6,3s2,3p6,3d10), is given by Fcore ðr;Se;28Þ ¼ F_ 1, 0 ðr;SeÞ þ F_ 2, 0 ðr;SeÞ þ F_ 2, 1 ðr;SeÞ þ F_ 3, 0 ðr;SeÞ þ F_ 3, 1 ðr;SeÞ þ F_ 3, 2 ðr;SeÞ
given by Fcore ðr;Se;18Þ ¼ F_ 1, 0 ðr;SeÞ þ F_ 2, 0 ðr;SeÞ þ F_ 2, 1 ðr;SeÞ þ FF 3, 0 ðr;SeÞ þ F_ 3, 1 ðr;SeÞ _
Because each subshell density for a given atom is expanded starting with the same set of primitive s-type Gaussian functions, one obtains a very compact, efficient expression for a total core density of the atom. The total charge density distribution of a molecule whose firstorder density matrix γ1(r,r0 ) was calculated using a Hamiltonian with ECPs can be simply approximated by adding the appropriate core electron density Fcore(r;ZA,Ncore,A) for each atom A whose Ncore,A electron core was modeled using an ECP to the electron density determined by γ1(r,r0 ): FðrÞ ¼ γ1 ðr;r0 Þjr ¼ r0 þ
ð6Þ while an 18 electron medium-core density Fcore(r;Se,18) for selenium atoms, corresponding to (1s2,2s2,2p6,3s2,3p6), is
ð7Þ
∑A Fcore ðr;ZA ;N core, A Þ
ð8Þ
Note that when using eq 8, it is also necessary to increment the 12883
dx.doi.org/10.1021/jp2040086 |J. Phys. Chem. A 2011, 115, 12879–12894
The Journal of Physical Chemistry A
ARTICLE
Table 3. QTAIM Atomic and Bonding Properties of ADSbO Using Different Methodsa property
A
B
C
D
E
Table 3. Continued property
F
q(Sb)
1.210
1.156
1.180
1.213
1.180
na
q(N) q(O)
1.264 1.124
1.240 1.108
1.242 1.117
1.261 1.128
1.242 1.117
na na
q(C1)
0.382
0.380
0.381
0.382
0.381
0.382
q(C2)
0.677
0.678
0.678
0.679
0.678
0.679
q(H1)
0.046
0.046
0.045
0.045
0.045
0.045
q(H2)
0.046
0.046
0.044
0.045
0.044
0.045
Vol(Mol)
1055.3
1046.5
1050.1
1053.3
1050.1
1053.3
Vol(Sb)
252.6
245.3
246.7
251.6
246.7
na
Vol(N) Vol(O)
88.3 114.0
87.6 113.5
88.0 114.3
88.3 113.6
88.0 114.3
na na
Vol(C1)
78.2
78.2
78.2
78.2
78.2
78.2
Vol(C2)
67.9
67.9
67.9
67.9
67.9
67.9
Vol(H1)
48.7
48.7
48.8
48.7
48.8
48.7
Vol(H2)
48.4
48.5
48.5
48.4
48.5
48.4
|μ(Mol)|
1.114
1.072
0.979
1.024
0.977
1.023
|μintra(Sb)|
1.584
1.454
1.411
1.514
1.411
na
|μintra(N)| |μintra(O)|
0.122 0.152
0.099 0.154
0.090 0.149
0.103 0.155
0.090 0.149
na na
|μintra(C1)|
0.655
0.653
0.654
0.654
0.654
0.654
|μintra(C2)|
0.754
0.751
0.753
0.754
0.753
0.754
|μintra(H1)|
0.138
0.138
0.139
0.138
0.139
0.138
|μintra(H2)|
0.136
0.136
0.136
0.135
0.136
0.135
Qxx(Sb)
9.248
8.986
9.096
9.368
9.096
na
Qxx(N)
1.921
1.876
1.882
1.954
1.882
na
Qxx(O) Qxx(C1)
0.946 3.669
0.936 3.672
0.951 3.674
0.962 3.675
0.951 3.674
na 3.674
Qxx(C2)
2.747
2.738
2.753
2.751
2.753
2.751
Qxx(H1)
0.277
0.276
0.278
0.279
0.278
0.279
Qxx(H2)
0.278
0.278
0.279
0.281
0.279
0.281
LI(Sb)
48.430
48.478
48.456
48.522
48.456
na
LI(Ni)
6.422
6.392
6.397
6.407
6.397
na
LI(O)
8.013
7.994
8.006
8.012
8.006
na
LI(C1) LI(C2)
3.687 3.449
3.688 3.448
3.686 3.447
3.686 3.447
3.686 3.447
3.686 3.447
LI(H1)
0.406
0.406
0.407
0.407
0.407
0.407
LI(H2)
0.415
0.415
0.416
0.415
0.416
0.416
DI(Sb|N)
0.852
0.862
0.863
0.830
0.863
na
DI(Sb|O)
0.665
0.665
0.663
0.640
0.663
na
DI(N|C1)
1.121
1.124
1.124
1.123
1.124
na
DI(O|C2)
1.115
1.121
1.118
1.117
1.118
na
DI(C1|C2) DI(C1|H1)
1.387 0.938
1.383 0.938
1.385 0.938
1.385 0.938
1.385 0.938
1.385 0.938
DI(C2|H2)
0.919
0.919
0.919
0.919
0.919
0.919
K(N)
55.223
55.418
55.205
55.215
55.205
na
K(O)
75.576
75.941
75.566
75.579
75.566
na
K(C1)
37.700
37.808
37.699
37.699
37.699
37.699
K(C2)
37.517
37.625
37.516
37.515
37.516
37.515
K(H1)
0.601
0.601
0.602
0.602
0.602
0.601
K(H2) Nleak(Sb)
0.605 na
0.605 na
0.605 0.0000
0.605 0.0569
0.605 0.0000
0.605 na
Fb(Sb|N)
0.1015
0.1023
0.1029
0.1052
0.1029
na
Fb(Sb|O)
0.0775
0.0780
0.0782
0.0799
0.0782
na
A
B
C
D
E
F
Fb(N|C1) Fb(O|C2)
0.3042 0.3286
0.3042 0.3289
0.3043 0.3287
0.3042 0.3286
0.3043 0.3287
0.3042 0.3286
Fb(C1|Cb)
0.3241
0.3241
0.3240
0.3242
0.3240
0.3242
Fb(C1|H1)
0.2817
0.2817
0.2817
0.2817
0.2817
0.2817
Fb(C2|H2)
0.2836
0.2835
0.2835
0.2835
0.2835
0.2835
r2Fb(Sb|N) 0.2614
0.2462
0.2648
0.2425
0.2648
na
r2Fb(Sb|O) 0.2512
0.2450
0.2580
0.2327
0.2580
na
r2Fb(N|C1) 0.6321 0.6331 0.6372 0.6351 0.6372 0.6351 r2Fb(O|C2) 0.2294 0.2426 0.2331 0.2326 0.2331 0.2326 r2Fb(C1|C2) 0.9284 0.9290 0.9288 0.9298 0.9288 0.9298 r2Fb(C1|H1) 0.9666 0.9660 0.9662 0.9663 0.9662 0.9663 r2Fb(C2|H2) 0.9874 0.9871 0.9870 0.9868 0.9870 0.9868 εb(Sb|N)
0.2928
0.2922
0.2826
0.2895
0.2826
na
εb(Sb|O)
0.1829
0.1837
0.1780
0.1830
0.1780
na
εb(N|C1)
0.1779
0.1800
0.1786
0.1785
0.1786
0.1785
εb(O|C2)
0.0326
0.0332
0.0336
0.0336
0.0336
0.0336
εb(C1|Cb) εb(C1|H1)
0.3483 0.0529
0.3455 0.0528
0.3471 0.0527
0.3465 0.0529
0.3471 0.0527
0.3465 0.0529
εb(C2|H2)
0.0392
0.0386
0.0393
0.0398
0.0393
0.0398
rb(Sb|N)
1.9654
1.9675
1.9731
1.9750
1.9731
na
rb(Sb|O)
2.0727
2.0773
2.0737
2.0742
2.0737
na
Rb(Sb|N)
3.9715
3.9715
3.9723
3.9715
3.9715
na
Rb(Sb|O)
4.0941
4.0941
4.0942
4.0948
4.0939
na
a
Geometry for all methods is the equilibrium geometry at B3LYP/ 6-311G(d,p)[H,C,N,O],cc-pVTZ-PP(S)[Sb], 28-electron ECP for Sb. A. B3LYP/6-311G(d,p)[H,C,N,O],DGDZVP+f[Sb], all-electron, nonrelativistic. B. B3LYP/6-311G(d,p)[H,C,N,O],DGDZVP+f[Sb], all-ellectron, relativistic (DKH2). C. B3LYP/6-311G(d,p)[H,C,N,O],cc-pVTZPP(S)[Sb], 28-electron ECP for Sb, 28-electron atomic core density on Sb. D. B3LYP/6-311G(d,p)[H,C,N,O],cc-pVTZ-PP(L)[Sb], 46-electron ECP for Sb, 46-electron atomic core density on Sb. E. B3LYP/6311G(d,p)[H,C,N,O],cc-pVTZ-PP(S)[Sb], 28-electron ECP for Sb, 28-electron tight core density on Sb. F. B3LYP/6-311G(d,p)[H,C,N,O], cc-pVTZ-PP(L)[Sb], 46-electron ECP for Sb, 46-electron tight core density on Sb.
nuclear charge of each atom A by Ncore,A, because in a calculation using ECPs effective nuclear charges ZA Ncore,A are used. As an alternative to the more realistic atomic core electron densities afforded by the above procedure, it is useful for comparison purposes to also consider representing a core electron density using a simplistic, tight core density function Fcoret(r;ZA,Ncore,A) such as the following: Fcore t ðr;ZA ;N core, A Þ ¼ 8N c exp½ 4πðr R A Þ2
ð9Þ
This simplistic single Gaussian function integrates to the number of ECP-modeled core electrons Ncore,A for atom A and is sufficient to restore the otherwise-absent local maximum of the electron density at the atom’s nucleus, but relative to Fc(r;ZA,Nc,A), it is tightly localized near the nucleus of the atom and, at least for third row atoms and higher, will not significantly “leak” density to the valence shell or to neighboring atoms.
’ CALCULATIONS AND SOFTWARE All ab initio calculations as well as all calculations of the fitted electron densities were done with a development version of Gaussian24 and with GaussView 5.0.9.38 The fitted atomic 12884
dx.doi.org/10.1021/jp2040086 |J. Phys. Chem. A 2011, 115, 12879–12894
The Journal of Physical Chemistry A
ARTICLE
Figure 2. Contour maps of the Laplacian of the electron density in the σh symmetry plane of ADSbO. Thin, solid, dark blue curves are positive Laplacian contours while thin, dashed dark red curves are negative Laplacian contours. Also shown are BCPs (small green spheres), RCPs (small red spheres), bond paths and interatomic surface paths. (a) All-electron nonrelativistic wave function. (b) All-electron relativistic wave function. (c) Wave function from a medium-core, 28-electron ECP on Sb calculation augmented with an Sb atomic small-core electron density. (d) Wave function from a large-core, 46-electron ECP on Sb calculation augmented with an Sb atomic large-core electron density.
subshell density data for all elements of the periodic table have been included in Gaussian.24 They are automatically combined as needed to produce appropriate total atomic core densities that are included in “extended wavefunction files”39 written by Gaussian when the wave function is based on a calculation for which at least one ECP was used. All QTAIM and other analyses of the extended wave functions produced by Gaussian, with and without additional core densities, were done with AIMAll.40 Atomic integration errors were less than 0.0001 au in most cases and less than 0.001 au in all cases, as measured both by the atomic Lagrangian values and by the additivity of atomic charges and kinetic energies. For consistency, all ab initio calculations in this work use the B3LYP DFT method2528 both for all-electron calculations and for calculations using ECPs. All-electron calculations were done both with and without scalar-relativistic
corrections for comparison with each other and with corresponding results from formally nonrelativistic calculations using ECPs. Scalar-relativistic calculations were done using the DouglasKrollHess second-order scalar-relativistic Hamiltonian (DKH2).3033
’ RESULTS The quantum theory of atoms in molecules (QTAIM)516 was used to analyze the electron density of several molecular systems containing heavy atoms to quantify local and regional differences in the electron density distributions from allelectron scalar-relativistic B3LYP and nonrelativistic B3LYP; from small-core, medium-core, and large-core ECPs with nonrelativistic B3LYP; and from different methods of representing (or not representing) the electron density contributions from the missing core electrons associated with ECPs. 12885
dx.doi.org/10.1021/jp2040086 |J. Phys. Chem. A 2011, 115, 12879–12894
The Journal of Physical Chemistry A
ARTICLE
QTAIM analyses provide properties for uniquely and physically quantifying local and regional properties of molecules.516 QTAIM properties are helpful not only for interpretation but also for prediction, when used as QSAR descriptors, for example. 4145 QTAIM properties can also be useful metrics for quantifying local and regional differences between the same molecule using different computational Table 4. QTAIM Atomic and Bonding Properties of Ferrocene (D5d) Using Different Methodsa property
A
B
C
D
E
F
q(Fe)
0.788
0.798
0.797
0.612
0.797
0.781
q(C)
0.113
0.114
0.117
0.094
0.117
0.110
q(H)
0.034
0.034
0.037
0.032
0.037
0.032
Vol(Mol) Vol(Fe)
1351.6 76.9
1351.5 77.3
1346.1 75.2
1344.1 72.0
1346.1 75.2
1343.9 66.2
Vol(C)
77.9
77.9
78.0
77.9
78.0
78.4
Vol(H)
49.6
49.5
49.1
49.4
49.1
49.4
|μintra(Fe)|
0
0
0
0
0
0
|μintra(C)|
0.173
0.173
0.190
0.180
0.190
0.168
|μintra(H)|
0.140
0.140
0.139
0.138
0.139
0.138
Qaa(Fe)
2.252
2.293
2.162
1.763
2.162
1.755
Qaa(C) Qaa(H)
2.825 0.257
2.827 0.256
2.869 0.252
2.796 0.253
2.869 0.252
2.903 0.253
Qbb(Fe)
2.251
2.293
2.162
1.763
2.162
1.755
Qbb(C)
1.070
1.073
1.071
1.013
1.071
1.086
Qbb(H)
0.194
0.194
0.198
0.200
0.198
0.200
Qcc(Fe)
4.503
4.586
4.324
3.526
4.324
3.510
Qcc(C)
1.755
1.754
1.798
1.783
1.798
1.816
Qcc(H)
0.451
0.450
0.450
0.453
0.450
0.453
LI(Fe) LI(C)
22.856 4.011
22.830 4.013
22.915 4.014
23.323 3.997
22.915 4.014
23.196 4.012
LI(H)
0.415
0.415
0.413
0.418
0.413
0.418
DI(Fe|C)
0.450
0.453
0.440
0.404
0.440
0.390
DI(C|C)
1.236
1.235
1.240
1.248
1.240
1.257
DI(C|H)
0.966
0.966
0.966
0.966
0.966
0.967
Nleak(Fe)
na
na
0.0000
0.0307
0.0000
0.0000
Fb(Fe|C)
0.0832
0.0836
0.0805
0.0795
0.0805
0.0727
Fb(C|C) Fb(C|H)
0.2884 0.2820
0.2883 0.2820
0.2883 0.2820
0.2887 0.2818
0.2883 0.2820
0.2887 0.2818
Fr3
0.0812
0.0815
0.0784
0.0774
0.0784
0.0712
Fr5
0.0496
0.0496
0.0514
0.0498
0.0514
0.0497
Fc
0.0493
0.0493
0.0498
0.0467
0.0498
0.0458
r2Fb(Fe|C) 0.2880
0.2844
0.3117
0.3160
0.3116
0.2397
r2Fb(C|C) 0.7188 0.7182 0.7151 0.7162 0.7151 0.7162 r Fb(C|H) 0.9727 0.9727 0.9728 0.9701 0.9728 0.9701 2
r2Fr3 r2Fr5
0.0812 0.0496
0.0815 0.0496
0.0784 0.0514
0.0774 0.0498
0.0784 0.0514
0.0712 0.0497
r2Fc
0.0492
0.0493
0.0498
0.0466
0.0498
0.0458
εb(Fe|C)
2.0851
2.0530
2.0586
2.3908
2.0586
2.4070
εb(C|C)
0.2563
0.2564
0.2554
0.2617
0.2554
0.2617
εb(C|H)
0.0228
0.0227
0.0219
0.0240
0.0219
0.0240
rb(Fe|C)
1.9463
1.9444
1.9623
1.9817
1.9623
1.9002
rb(C|Fe)
1.9693
1.9712
1.9531
1.9322
1.9531
2.0180
rb(C|C) rb(C|H)
1.3471 1.3074
1.3471 1.3076
1.3471 1.3086
1.3470 1.3058
1.3471 1.3086
1.3470 1.3058
Table 4. Continued property
A
B
C
D
E
F
rb(H|C) Rb(Fe|C)
0.7327 3.9496
0.7324 3.9496
0.7315 3.9459
0.7342 3.9468
0.7315 3.9355
0.7342 3.9945
Rb(C|C)
2.6952
2.6952
2.6952
2.6946
2.6952
2.6946
Rb(C|H)
2.0121
2.0121
2.0121
2.0122
2.0121
2.0122
a
Geometry for all methods is the equilibrium geometry at B3LYP/ 6-311G(d,p)[H,C,Fe], all-electron, nonrelativistic. A. B3LYP/6-311G (d,p)[H,C,Fe], all-electron, nonrelativistic. B. B3LYP/6-311G(d,p)[H,C,Fe], all-electron, relativistic (DKH2). C. B3LYP/6-311G(d,p)[H,C],LANL2DZ[Fe], 10-electron ECP for Fe, 10-electron atomic core density on Fe. D. B3LYP/6-311G(d,p)[H,C],LANL1DZ[Fe], 18-electron ECP for Fe, 18-electron atomic core density on Fe. E. B3LYP/6-311G(d,p)[H,C],LANL2DZ[Fe], 10-electron ECP for Fe, 10-electron tight core density on Fe. F. B3LYP/6-311G(d,p)[H,C],LANL1DZ[Fe], 18-electron ECP for Fe, 18-electron tight core density on Fe.
Figure 3. Molecular graph of ferrocene molecule calculated from a wave function of a B3LYP calculation using a 6-311G(d,p) basis for C and H and the LANL2DZ ECP and basis for Fe. The total electron density from the wave function is augmented with an Fe small-core, 10-electron atomic core density. Small green spheres are bond critical points. Small red spheres are ring critical points. Small blue spheres are cage critical points. Dark red paths are unique gradient paths of the electron density connecting RCPs to associated BCPs.
methods, which is the primary use to which they shall be put to here. QTAIM has been amply applied and reviewed in the literature,516 but for clarity the QTAIM property symbols used in this paper are briefly summarized in Table 1, most of which are identical to those given in Bader’s 1990 QTAIM book.5 GermaniumAdamantane (Ge10H16). An example that clearly illustrates the problem of missing core electron densities from ECP-based wave functions is Geadamantane (adamantane with all carbon atoms replaced by germanium atoms), which has Td symmetry. Figure 1 shows the QTAIM molecular graph and the symmetrically unique atomic charges and electron density values at bond critical points 12886
dx.doi.org/10.1021/jp2040086 |J. Phys. Chem. A 2011, 115, 12879–12894
The Journal of Physical Chemistry A
ARTICLE
atomic core density functions, respectively. The geometry in all cases is the all-electron B3LYP/6-311G(d,p) nonrelativistic equilibrium geometry. In the small-core ECP case with no core density augmentation, Figure 1c, one sees numerous spurious, bond, ring, cage, and non-nuclear attractor critical points of the electron density, but only in the inner core region of the Ge atoms, near the nuclei. These spurious critical points might cause difficulties for some QTAIM topological analysis and atomic integration software, but because they are not topologically connected to genuine critical points outside the inner core (in this case), the problems are not Table 5. Difference between QTAIM Atomic and Bonding Properties of Transplatin and Cisplatin Using Different Methodsa property
Figure 4. Molecular graphs of the platinum compounds cisplatin and transplatin together with platinum interatomic surface paths in the plane of the non-hydrogen atoms. Electron density isosurfaces (0.001) colormapped with the total electrostatic potential are also shown. The electron densities were calculated from ECP-based B3LYP/6-311++G(2d,2p)[H,N,Cl],LANL08(f)[Pt] equilibrium geometry wave functions augmented with a 60-electron Pt atomic core density.
(BCPs) of this molecule as determined from the electron density distribution obtained using five different methods. Electron density contours in a σd symmetry plane are also shown. For Figure 1a, the electron density distribution is from an allelectron B3LYP/6-311G(d,p)4650 nonrelativistic wave function. In Figure 1b, the electron density is from an all-electron B3LYP/6-311G(d,p) relativistic wave function. In Figure 1c, the electron density is from a B3LYP/6-311G(d,p)[C,H], cc-pVTZPP[Ge]51 wave function in which small-core, 10-electron ECPs were used for the Ge atoms. In Figure 1d the electron density was calculated from a B3LYP/6-311G(d,p)[C,H], SDB-cc-pVTZ [Ge]52 wave function in which large-core, 28-electron ECPs were used for the Ge atoms. Parts e and f of Figure 1 are the QTAIM results obtained when the wave functions corresponding to Figure 1c,d are augmented with 10-electron and 28-electron Ge
A
B
C
D
q(Pt)
0.017
0.014
0.013
0.013
q(N)
0.015
0.010
0.012
0.012
q(Hi)
0.018
0.017
0.018
0.018
q(Ho)
0.022
0.021
0.021
0.021
q(Cl)
0.050
0.043
0.043
0.043
|μ(Mol)|
4.240
4.205
4.220
4.220
|μintra(Pt)|
0.055
0.082
0.075
0.075
|μintra(N)| |μintra(Hi)|
0.017 0.001
0.012 0.002
0.014 0.001
0.013 0.001
|μintra(Ho)|
0.007
0.006
0.006
0.006
|μintra(Cl)|
0.076
0.079
0.095
0.095
E(Mol)
0.018
0.021
0.021
0.021
K(Mol)
0.017
0.029
0.007
0.007
K(Pt)
0.033
0.045
0.009
0.009
K(N)
0.006
0.007
0.007
0.007
K(Hi) K(Ho)
0.013 0.011
0.013 0.011
0.013 0.011
0.013 0.011
K(Cl)
0.005
0.007
0.006
0.007
LI(Pt)
0.027
0.020
0.018
0.018
LI(N)
0.030
0.021
0.024
0.023
LI(Hi)
0.014
0.014
0.014
0.014
LI(Ho)
0.012
0.012
0.012
0.012
LI(Cl)
0.098
0.083
0.085
0.085
DI(Pt|N) DI(Pt|Cl)
0.041 0.033
0.037 0.031
0.038 0.033
0.038 0.033
DI(N|Hi)
0.040
0.040
0.041
0.041
DI(N|Ho)
0.024
0.022
0.022
0.022
Nleak(Pt)
na
na
0.0000
0.0000
Fb(Pt|N)
0.0113
0.0109
0.0104
0.0104
Fb(Pt|Cl)
0.0057
0.0058
0.0057
0.0057
Fb(N|Hi)
0.0057
0.0057
0.0058
0.0058
Fb(N|Ho) rb(Pt|N)
0.0003 0.0439
0.0004 0.0470
0.0004 0.0479
0.0004 0.0481
rb(N|Pt)
0.0348
0.0319
0.0310
0.0309
rb(Pt|Cl)
0.0190
0.0224
0.0244
0.0245
rb(Cl|Pt)
0.0253
0.0217
0.0198
0.0198
rb(N|Hi)
0.0175
0.0175
0.0178
0.0178
rb(Hi|N)
0.0075
0.0074
0.0077
0.0077
rb(N|Ho)
0.0090
0.0086
0.0086
0.0086
rb(Ho|N) Kb(Pt|N)
0.0074 0.0072
0.0070 0.0075
0.0071 0.0060
0.0071 0.0060
12887
dx.doi.org/10.1021/jp2040086 |J. Phys. Chem. A 2011, 115, 12879–12894
The Journal of Physical Chemistry A
ARTICLE
Table 5. Continued property
A
B
C
D
Kb(Pt|Cl) Kb(N|Hi)
0.0040 0.0003
0.0045 0.0002
0.0033 0.0003
0.0033 0.0003
Kb(N|Ho)
0.0089
0.0082
0.0081
0.0081
Gb(Pt|N)
0.0148
0.0162
0.0128
0.0130 0.0017
Gb(Pt|Cl)
0.0021
0.0025
0.0018
Gb(N|Hi)
0.0007
0.0006
0.0005
0.0005
Gb(N|Ho)
0.0026
0.0024
0.0023
0.0023
a
Geometry for all methods is the equilibrium geometry at B3LYP/6-311 ++G(2d,2p)[H,N,Cl],LANL08(f)[Pt], 60-electron ECP for Pt. A. B3LYP/6-311++G(2d,2p)[H,N,Cl],SARC-DKH[Pt] all-electron, nonrelativistic. B. B3LYP/6-311++G(2d,2p)[H,N,Cl],SARC-DKH[Pt], allelectron, relativistic (DKH2). C. B3LYP/6-311++G(2d,2p)[H,N,Cl], LANL08(f)[Pt], 60-electron ECP for Pt, 60-electron atomic core density on Pt. D. B3LYP/6-311++G(2d,2p)[H,N,Cl],LANL08(f)[Pt], 60-electron ECP for Pt, 60-electron tight core density on Pt.
fatal and can be avoided for most purposes by simply ignoring them. Comparison of the atomic charges and BCP electron density values of Figure 1b shows quite reasonable agreement with Figure 1a,b, the values lying between the nonrelativistic and relativistic all-electron values. In the large-core ECP case with no core density augmentation, Figure 1d, the definition of all atoms, including the hydrogens, and their connectivity has been qualitatively altered due to the presence of spurious non-nuclear attractor critical points in the valence regions. QTAIM analysis is meaningless in this case. For this reason, atomic charges and BCP electron density values cannot be shown in Figure 1c. Parts e and f of Figure 1 show that atomic core density augmentation of both the small-core and large-core ECP-based wave functions results in the correct topology of the total electron density and reasonable agreement of the atomic charges with the all electron cases. Not surprisingly, the atomic charge and BCP electron density data from the small-core case, Figure 1e, is in better agreement with the all-electron results than the data from the large-core case and they are also identical to the corresponding result obtained from no core density augmentation, Figure 1c. That the atomic charges and BCP electron density values in the small-core case are not closer to the all-electron results is due to the use of a formally nonrelativistic ECP Hamiltonian instead of an all-electron relativistic or nonrelativistic Hamiltonian and also the use of a somewhat different valence basis set. In the largecore case, possible “contaminating” effects of the augmentative core densities on the electron density in the valence region is also a factor. One measure of this is the quantity Nleak(A): Z ð10Þ N leak ðAÞ ¼ N core, A dr Fcore ðr;ZA ;N core, A Þ A
One expects this quantity to be zero when augmenting with small-core, medium-core and fully appropriate large-core densities. In the Geadamantane case, Nleak(Ge) is effectively zero with the small-core density augmentation but for large-core density augmentation is about 0.032 for the methylene-like Ge atoms and about 0.020 for the methine-like Ge atoms, which is consistent with the slightly higher positive charge of the Ge atoms in the large-core case compared to the small-core case and the all-electron cases.
Figure 5. Molecular graphs of Td (a) and C3v (b) isomers of tetrairidium dodecacarbonyl. The molecular graphs were derived from electron densities calculated from ECP-based B3LYP/cc-pVDZ[C,O], cc-pVDZ-PP[Ir] equilibrium geometry wave functions augmented with 60-electron Ir atomic core densities.
These results for Geadamantane are consistent with the results of Tiana et al.21 in their HartreeFock study of the effects of core density augmentation on QTAIM properties of GeH4. Titanium Dioxide. Table 2 shows some QTAIM atomic and bonding properties of titanium dioxide using four variants of the B3LYP method: nonrelativistic, all-electron B3LYP using the cc-pVDZ53,54 basis set for all atoms; relativistic, all-electron B3LYP using the cc-pVDZ basis for all atoms; B3LYP using the cc-PVDZ basis for the oxygen atoms and the LANL2DZ55 small-core, 10electron ECP and corresponding basis for the titanium atom; B3LYP using the cc-PVDZ basis for the oxygen atoms and the LANL1DZ56 large-core, 18-electron ECP and corresponding basis for the titanium atom. For the two methods employing an ECP for the titanium atom, data for three different methods of representing the core electron density of the titanium atom are shown: no core density, an atomic core density constructed from 12888
dx.doi.org/10.1021/jp2040086 |J. Phys. Chem. A 2011, 115, 12879–12894
The Journal of Physical Chemistry A
ARTICLE
Table 6. QTAIM Atomic and Bonding Properties of Tetrairidium Dodecacarbonyl (Td) Using Different Methodsa,b property
A
B
C
Table 7. QTAIM Atomic and Bonding Properties of Tetrairidium Dodecacarbonyl (C3v) Using Different Methodsa,b
D
property
A
B
C
D
q(Ir)
0.345
0.383
0.340
0.350
q(Ir1)
0.314
0.334
0.290
0.299
q(C) q(O)
1.052 1.167
1.041 1.168
1.056 1.169
1.052 1.169
q(Ir2) q(C7)
0.409 1.060
0.457 1.045
0.416 1.060
0.425 1.057
|μ(Mol)|
0
0
0
0
q(C5)
0.927
0.917
0.924
0.922
|μintra(Ir)|
0.478
0.5477
0.503
0.512
q(C17)
1.075
1.061
1.078
1.075
|μintra(C)|
1.355
1.395
1.373
1.378
q(C13)
1.075
1.054
1.070
1.067
|μintra(O)|
0.890
0.886
0.885
0.885
q(O8)
1.165
1.165
1.166
1.166
E(Mol)
68461.053
72598.021
1778.846
1778.846
q(O6)
1.163
1.158
1.154
1.154
K(Mol)
68210.391
102838.224
1492.345
1492.345
q(O18)
1.163
1.161
1.162
1.162
K(Ir) K(C)
16714.351 37.056
25371.327 37.041
34.902 37.036
34.899 37.037
q(O14) |μ(Mol)|
1.160 0.197
1.161 0.177
1.162 0.205
1.162 0.205
K(O)
75.693
75.701
75.692
75.692
|μintra(Ir1)|
0.459
0.572
0.532
0.542
DI(Ir|Ir)
0.620
0.695
0.694
0.693
|μintra(Ir2)|
0.241
0.289
0.242
0.250
DI(Ir|C)
1.127
1.218
1.217
1.217
|μintra(C7)|
1.362
1.401
1.380
1.385
DI(C|O)
1.549
1.531
1.534
1.534
|μintra(C5)|
1.338
1.338
1.323
1.326
Nleak(Ir)
na
na
0.0000
0.0000
|μintra(C17)|
1.326
1.368
1.342
1.347
Fb(Ir|Ir)
0.0610
0.0672
0.0659
0.0658
|μintra(C13)|
1.352
1.393
1.369
1.374
Fb(Ir|C) Fb(C|O)
0.1562 0.4551
0.1582 0.4550
0.1598 0.4549
0.1591 0.4549
|μintra(O8)| |μintra(O6)|
0.896 0.789
0.891 0.774
0.890 0.770
0.890 0.770
Fr
0.0436
0.0465
0.0464
0.0464
|μintra(O18)|
0.890
0.884
0.884
0.884
Fc
0.0410
0.0427
0.0425
0.0425
|μintra(O14)|
0.898
0.891
0.890
0.890
V b (Ir|Ir)
0.0552
0.0548
0.0504
0.0504
E(Mol)
68461.100
72598.020
1778.842
1778.842
V b (Ir|C)
0.2687
0.2717
0.2633
0.2655
K(Mol)
68210.528
102838.054
1492.629
1492.629
V b (C|O)
1.8679
1.8691
1.8676
1.8676
K(Ir1)
16714.363
25371.320
34.949
34.946
Vr
0.0320
0.0341
0.0341
0.0341
K(Ir2)
16714.384
25371.293
34.986
34.982
Vc
0.0267
0.0297
0.0297
0.0297
K(C7) K(C5)
37.054 37.101
37.037 37.066
37.032 37.073
37.033 37.074
K(C17)
37.062
37.043
37.036
37.038
K(C13)
37.059
37.045
37.042
37.043
K(O8)
75.694
75.702
75.693
75.693
K(O6)
75.656
75.666
75.655
75.655
K(O18)
75.689
75.695
75.686
75.686
K(O14)
75.690
75.699
75.690
75.690
DI(Ir1|Ir2) DI(Ir1|C7)
0.585 1.112
0.661 1.204
0.660 1.203
0.660 1.203
a
Geometry for all methods is the equilibrium geometry at B3LYP/ccpVDZ[C,O],cc-pVDZ-PP[Ir], 60-electron ECP for Ir. A. B3LYP/ccpVDZ[C,O],SARC-DKH[Ir] all-electron, nonrelativistic. B. B3LYP/ccpVDZ[C,O],SARC-DKH[Ir] all-electron, relativistic (DKH2). C. B3LYP/cc-pVDZ[C,O],cc-pVDZ-PP[Ir], 60-electron ECP for Ir, 60electron atomic core density on Ir. D. B3LYP/cc-pVDZ[C,O],ccpVDZ-PP[Ir], 60-electron ECP for Ir, 60-electron tight core density on Ir. b See Supporting Information Table S5 for a comparison of more properties.
fitted core density subshells of a free titanium atom, and the tight core density function defined by eq 9. The geometry in all cases is the all-electron, nonrelativistic B3LYP/cc-pVDZ equilibrium geometry. The data in Table 2 show that there is just a small difference between the all-electron nonrelativistic and relativistic results, which is not too surprising because titanium is only an early firstrow transition metal. The small-core LANL2DZ results are in good agreement with the all-electron results. The use of no core density or a 10-electron tight core density, eq 9, for the Ti small core is essentially identical to the use of a 10-electron Ti atomic core density. Interestingly, the LANL2DZ QTAIM results are generally closer to the nonrelativistic all-electron QTAIM results than the corresponding relativistic ones. This may be related to the fact that the LANL2DZ ECPs for the first-row transition metals are based on nonrelativistic HartreeFock calculations. In contrast, the large-core LANL1DZ results are quantitatively correct for some properties but are only qualitatively correct for others, compared to the all-electron results. When an 18-electron atomic core density is used for the large core of Ti, the net
DI(Ir2|Ir3)
0.398
0.457
0.453
0.452
DI(Ir2|C5)
0.689
0.740
0.738
0.738
DI(Ir2|C17)
1.194
1.287
1.285
1.285
DI(Ir2|C13)
1.153
1.244
1.241
1.241
DI(C7|O8)
1.559
1.539
1.543
1.543
DI(C5|O6)
1.502
1.486
1.493
1.493
DI(C17|O18) DI(C13|O14)
1.552 1.563
1.531 1.543
1.536 1.547
1.536 1.547
Nleak(Ir1)
na
na
0.0000
0.0000
Nleak(Ir2)
na
na
0.0000
0.0000
Fb(Ir1|Ir2)
0.0584
0.0650
0.0637
0.0636
Fb(Ir1|C7)
0.1532
0.1553
0.1568
0.1561
Fb(Ir2|Ir3)
0.0502
0.0546
0.0537
0.0536
Fb(Ir2|C5)
0.1084
0.1123
0.1113
0.1111
Fb(Ir2|C17) Fb(Ir2|C13)
0.1643 0.1552
0.1666 0.1569
0.1687 0.1585
0.1679 0.1577
Fb(C7|O8)
0.4562
0.4560
0.4559
0.4559
Fb(C5|O6)
0.4376
0.4372
0.4367
0.4367
12889
dx.doi.org/10.1021/jp2040086 |J. Phys. Chem. A 2011, 115, 12879–12894
The Journal of Physical Chemistry A
ARTICLE
Table 7. Continued property
A
B
C
D
Fb (C17|O18) Fb(C13|O14)
0.4560 0.4567
0.4558 0.4565
0.4558 0.4565
0.4558 0.4565
Fr(Ir1|Ir2|Ir3)
0.0765
0.0792
0.0713
0.0712
Fr(Ir2|Ir3|Ir4)
0.0501
0.0454
0.0430
0.0430
Fr(Ir2|Ir4|C5)
0.1400
0.1290
0.1295
0.1287
Fc
0.0376
0.0407
0.0404
0.0404
V b (Ir1|Ir2)
0.0503
0.0506
0.0464
0.0464
V b (Ir1|C7)
0.2639
0.2659
0.2578
0.2599
V b (Ir2|Ir3) V b (Ir2|C5)
0.0489 0.1432
0.0484 0.1419
0.0447 0.1334
0.0447 0.1338
V b (Ir2|C17)
0.2892
0.2906
0.2829
0.2854
V b (Ir2|C13)
1.2002
1.2080
1.2097
1.2095
V b (C7|O8)
1.8741
1.8761
1.8743
1.8743
V b (C5|O6)
1.7534
1.7576
1.7529
1.7529
V b (C17|O18)
1.8739
1.8764
1.8749
1.8749
V b (C13|O14)
1.8785
1.8806
1.8786
1.8786
V r (Ir1|Ir2|Ir3) V r (Ir2|Ir3|Ir4)
0.0281 0.0244
0.0296 0.0249
0.0310 0.0257
0.0310 0.0257
V r (Ir2|Ir4|C5)
0.0499
0.0500
0.0461
0.0461
Vc
0.0235
0.0250
0.0268
0.0269
a
Geometry for all methods is the equilibrium geometry at B3LYP/ccpVDZ[C,O],cc-pVDZ-PP[Ir], 60-electron ECP for Ir. A. B3LYP/ccpVDZ[C,O],SARC-DKH[Ir] all-electron, nonrelativistic. B. B3LYP/ccpVDZ[C,O],SARC-DKH[Ir] all-electron, relativistic (DKH2). C. B3LYP/cc-pVDZ[C,O],cc-pVDZ-PP[Ir], 60-electron ECP for Ir, 60electron atomic core density on Ir. D. B3LYP/cc-pVDZ[C,O],cc-pVDZPP[Ir], 60-electron ECP for Ir, 60-electron tight core density on Ir. b See Supporting Information Table S6 for a comparison of more properties.
electron population of Ti is larger by more than 0.3 electrons and the volume of the Ti atom is about 9% larger than the all-electron case. When an 18-electron tight core density is used for the largecore Ti atom, the volume of the Ti atom is about 25% too small and the electron population is correspondingly too low by about 0.15 electrons. The number of “leaked” core electrons from the titanium to the oxygen atoms, Nleak(Ti), is about 0.123 in the large-core LANL1DZ case when the core density is represented using an atomic core density, but Nleak(Ti) is zero in all other cases. The distance from a Ti|O bond critical point to the Ti nucleus is too large by about 0.06 au when using LANL1DZ with the atomic core density but is about 0.20 au too small when using the tight core density, which is consistent with many of the other differences in atomic and bonding properties. The intra-atomic dipole moment magnitude of the Ti atom, |μintra(Ti)|, is much smaller in the large-core case with an atomic core density, which might suggest that the dipole polarization contribution from the outer core shell of the Ti atom is significant. Because of the design of the tight core density, the bond critical point properties for both the small-core and large-core cases when no core density is present are identical to those of the corresponding tight core density cases, as would the atomic properties be. In this case (but not necessarily in general) the atomic properties in the no core density cases can actually determined by simply ignoring all of the spurious critical points, which are all in the core of the Ti atom and not connected to the critical points in the valence region. As was seen with the large-core case of Ge10H16 and as will be seen in later examples, TiO2 is not necessarily typical and
Table 8. Difference between QTAIM Atomic and Bonding Properties of Tetrairidium Dodecacarbonyl (C3v) and Tetrairidium Dodecacarbonyl (Td) Using Different Methodsa,b property
A
B
C
D
q(Ir1)
0.031
0.049
0.050
0.051
q(Ir2)
0.064
0.074
0.076
0.075
q(C7)
0.008
0.004
0.004
0.005
q(C5)
0.125
0.124
0.132
0.130
q(C17) q(C13)
0.023 0.023
0.020 0.013
0.022 0.014
0.023 0.015
q(O8)
0.002
0.003
0.003
0.003
q(O6)
0.004
0.010
0.015
0.015
q(O18)
0.004
0.007
0.007
0.007
q(O14)
0.007
0.007
0.007
0.007
|μ(Mol)|
0.197
0.177
0.205
0.205
|μintra(Ir1)|
0.019
0.025
0.029
0.030
|μintra(Ir2)| |μintra(C7)|
0.237 0.007
0.259 0.006
0.261 0.007
0.262 0.007
|μintra(C5)|
0.017
0.057
0.050
0.052
|μintra(C17)|
0.029
0.027
0.031
0.031
|μintra(C13)|
0.003
0.002
0.004
0.004
|μintra(O8)|
0.006
0.005
0.005
0.005
|μintra(O6)|
0.101
0.112
0.115
0.115
|μintra(O18)|
0.000
0.002
0.001
0.001
|μintra(O14)| E(Mol)
0.008 0.047
0.005 0.001
0.005 0.004
0.005 0.004
K(Mol)
0.137
0.170
0.284
0.284
K(Ir1)
0.012
0.007
0.047
0.047
K(Ir2)
0.033
0.034
0.084
0.083
K(7)
0.002
0.004
0.004
0.004
K(C5)
0.045
0.025
0.037
0.037
K(C17)
0.006
0.002
0.000
0.001
K(C13) K(O8)
0.003 0.001
0.004 0.001
0.006 0.001
0.006 0.001
K(O6)
0.037
0.035
0.037
0.037
K(O18)
0.004
0.006
0.006
0.006
K(O14)
0.003
0.002
0.002
0.002
DI(Ir1|Ir2)
0.035
0.034
0.034
0.033
DI(Ir1|C7)
0.015
0.014
0.014
0.014
DI(Ir2|Ir3)
0.222
0.238
0.241
0.241
DI(Ir2|C5) DI(Ir2|C17)
0.438 0.067
0.478 0.069
0.479 0.068
0.479 0.068
DI(Ir2|C13)
0.026
0.026
0.024
0.024
DI(C7|O8)
0.010
0.008
0.009
0.009
DI(C5|O6)
0.047
0.045
0.041
0.041
DI(C17|O18)
0.003
0.000
0.002
0.002
DI(C13|O14)
0.014
0.012
0.013
0.013
Nleak(Ir1)
na
na
0.0000
0.0000
Nleak(Ir2) Fb(Ir1|Ir2)
na 0.0026
na 0.0022
0.0000 0.0022
0.0000 0.0022
Fb(Ir1|C7)
0.0030
0.0029
0.0030
0.0030
Fb(Ir2|Ir3)
0.0108
0.0126
0.0122
0.0122
Fb(Ir2|C5)
0.0478
0.0459
0.0485
0.0480
Fb(Ir2|C17)
0.0081
0.0084
0.0089
0.0088
Fb(Ir2|C13)
0.0010
0.0013
0.0013
0.0014
Fb(C7|O8)
0.0011
0.0010
0.0010
0.0010
12890
dx.doi.org/10.1021/jp2040086 |J. Phys. Chem. A 2011, 115, 12879–12894
The Journal of Physical Chemistry A
ARTICLE
Table 8. Continued property
A
B
C
D
Fb(C5|O6) Fb (C17|O18)
0.0175 0.0009
0.0178 0.0008
0.0182 0.0009
0.0182 0.0009
Fb(C13|O14)
0.0016
0.0015
0.0016
0.0016
Fr(Ir1|Ir2|Ir3)
0.0329
0.0327
0.0249
0.0248
Fr(Ir2|Ir3|Ir4)
0.0065
0.0011
0.0034
0.0034
Fr(Ir2|Ir4|C5)
na
na
na
na
Fc
0.0034
0.0020
0.0021
0.0021
V b (Ir1|Ir2)
0.0050
0.0042
0.0040
0.0040
V b (Ir1|C7) V b (Ir2|Ir3)
0.0048 0.0063
0.0058 0.0064
0.0055 0.0057
0.0056 0.0057
V b (Ir2|C5)
0.1255
0.1298
0.1299
0.1317
V b (Ir2|C17)
0.0205
0.0189
0.0196
0.0199
V b (Ir2|C13)
0.9315
0.9363
0.9464
0.9440
V b (C7|O8)
0.0062
0.0070
0.0067
0.0067
V b (C5|O6)
0.1145
0.1115
0.1148
0.1148
V b (C17|O18)
0.0060
0.0073
0.0073
0.0073
V b (C13|O14) V r (Ir1|Ir2|Ir3)
0.0106 0.0039
0.0115 0.0045
0.0110 0.0031
0.0110 0.0031
V r (Ir2|Ir3|Ir4)
0.0076
0.0092
0.0084
0.0084
V r (Ir2|Ir4|C5)
na
na
na
na
Vc
0.0032
0.0047
0.0029
0.0029
a
Geometry for all methods is the equilibrium geometry at B3LYP/ccpVDZ[C,O],cc-pVDZ-PP[Ir], 60-electron ECP for Ir. A. B3LYP/ccpVDZ[C,O],SARC-DKH[Ir] all-electron, nonrelativistic. B. B3LYP/ccpVDZ[C,O],SARC-DKH[Ir] all-electron, relativistic (DKH2). C. B3LYP/cc-pVDZ[C,O],cc-pVDZ-PP[Ir], 60-electron ECP for Ir, 60electron atomic core density on Ir. D. B3LYP/cc-pVDZ[C,O],ccpVDZ-PP[Ir], 60-electron ECP for Ir, 60-electron tight core density on Ir. b See Supporting Information Table S7 for a comparison of more properties.
it can and does happen that the lack of a core density for an atom can lead to spurious critical points even outside of the core, which makes QTAIM analyses impossible without restoring the core density. Also noteworthy in Table 2 is that the molecular dipole moments show less variation among the different methods than the intra-atomic dipole moments. This is illustrative of the general point that atomic properties are generally more sensitive indicators of differences between molecules or between methods than are molecular properties. Similarly, local properties are generally more sensitive indicators of differences than atomic properties. ADSbO. Table 3 shows some QTAIM atomic and bonding properties of an ADSbO57 molecule containing the heavy main group element antimony using four variants of the B3LYP method: nonrelativistic, all-electron B3LYP using the 6-311G(d,p)4650 basis for the H, C, N, and O atoms and the DGDZVP58,59 basis with an additional shell of F basis functions for the antimony atom; relativistic, all-electron B3LYP using the 6-311G(d,p) basis for the H, C, N, and O atoms and the DGDZVP basis with an additional shell of F basis functions for the antimony atom; B3LYP using the 6-311G(d,p) basis for the H, C, N, and O atoms and the cc-pVTZ-PP(S)60 medium-core, 28-electron ECP and corresponding basis for the antimony atom; B3LYP using the 6-311G(d,p) basis for the H, C, N, and O atoms and the SDB-ccpVTZ52 large-core, 46-electron ECP and corresponding basis for the antimony atom. For the two B3LYP methods employing an
ECP for the antimony atom, data for two different methods of representing the core electron density of the antimony atom are shown: an Sb atomic core density constructed from the core density subshells of a free Sb atom; and the tight core density function defined by eq 9. The geometry in all cases is the equilibrium geometry of the medium-core ECP-based method, B3LYP/6-311G(d,p)[H,C,N,O],cc-pVTZ-PP(S)[Sb]. Both the medium-core and large-core cases show good agreement with the all-electron cases when augmenting with atomic core densities, the medium-core results giving a little better agreement, not surprisingly. The use of a tight core density with the medium-core ECP-based wave function gives results virtually identical to those obtained when augmenting with an atomic medium-core density. Importantly, the large-core case is an example where the use of a simplistic tight core density function, eq 9, is inadequate, and a more realistic atomic core density is necessary to obtain the correct topology and QTAIM properties. In contrast with the TiO2 results, the augmented small-core QTAIM results are generally closer to the relativistic all-electron QTAIM results than the corresponding nonrelativistic ones. This is a reflection of the greater importance of relativistic effects for the much higher nuclear charge of the Sb atom and the fact that the ECP for the Sb atom was based on relativistic calculations. Figure 2 compares maps of the Laplacian of the electron density together with critical points, bond paths, and interatomic surface paths in the σh symmetry plane for the all-electron wave functions and the ECP-based wave functions augmented with atomic core densities. Clearly the agreement between the maps reflects what the data in Table 3 show. Ferrocene (D5d). Table 4 shows some QTAIM atomic and bonding properties of a D5d ferrocene molecule (Figure 3) using four variants of the B3LYP DFT method: nonrelativistic, allelectron B3LYP using the 6-311G(d,p)4650 basis for all atoms; relativistic, all-electron B3LYP using the 6-311G(d,p) basis for all atoms; B3LYP using the 6-311G(d,p) basis for the H and C atoms and the LANL2DZ55 small-core, 10-electron ECP and corresponding basis for the iron atom; B3LYP using the 6-311G(d,p) basis for the H and C atoms and the LANL1DZ56 large-core, 18-electron ECP and corresponding basis for the iron atom. For the two B3LYP methods employing an ECP for the iron atom, data for two different methods of representing the core electron density of the iron atom are shown: an atomic core density constructed from the core density subshells of a free iron atom; and the tight core density function defined by eq 9. The geometry in all cases is the B3LYP/6-311G(d,p) all-electron, nonrelativistic equilibrium geometry. The data in Table 4 lead to conclusions similar to that for TiO2. The nonrelativistic and relativistic all-electron results are virtually identical. The small-core ECP results with either an atomic core density or tight core density for the Fe atom are essentially identical and in quite reasonable agreement with the all-electron results. The large-core ECP results are in qualitative agreement with the all-electron results but there is a noticeable difference between the atomic core density augmentation and the tight core density augmentation, with the former generally giving better agreement with the allelectron results. Transplatin and Cisplatin. Supporting Information Tables S4 and S4 compare QTAIM atomic and bonding properties for the platinum compound transplatin61 (C2h symmetry) and its anticancer isomer cisplatin61 (C2v symmetry) calculated 12891
dx.doi.org/10.1021/jp2040086 |J. Phys. Chem. A 2011, 115, 12879–12894
The Journal of Physical Chemistry A using wave functions from nonrelativistic and relativistic allelectron B3LYP calculations as well as from a B3LYP calculation using a 60-electron ECP for the platinum atom. In the latter case, augmentation with a Pt 60-electron atomic core density is compared to augmentation with a tight core density function, eq 9. The basis set used for the non-Pt atoms was 6-311++G(2d,2p).4650,62 The all-electron basis used for the Pt atom was SARC-DKH.63 The ECP basis and potential used for Pt was LANL08(f).55,56,6466 Figure 4 shows QTAIM molecular graphs, 0.001 electron density isosurfaces mapped with the total electrostatic potential and platinum interatomic surface paths in the plane of the nonhydrogen atoms for transplatin and cisplatin using the ECPbased wave functions augmented with Pt atomic core densities. Note that the ESP-mapped isodensity surfaces are shown for effect but the representation of the Pt core makes essentially no contribution to the determination of a 0.001 electron density isosurface or to the electrostatic potential on it because the neutral Pt core is confined well inside it. Table 5 shows the differences in QTAIM properties between transplatin and cisplatin to see how well differences in QTAIM properties between the isomers compare using different different wave functions and methods of representing the platinum core density for the ECP-based wave function. The zero values for Nleak(Pt) and the essentially exact agreement between the use of an atomic Pt core density versues a tight Pt core density ensures that the differences in QTAIM results between the ECP-based wave function augmented with an atomic Pt core density and the all-electron wave functions is due only to the wave functions. With the exception of the non-hydrogen atomic kinetic energies, the agreement between the absolute ECP-based QTAIM results and all-electron results is quite good for each molecule and the agreement between methods is even better when comparing differences between the isomers. That the nitrogen and chlorine absolute atomic kinetic energies are significantly different between the ECP-based methods and allelectron methods is partly due to basis set differences. However, the relative kinetic energies as well as all other properties of all nonplatinum atoms are in excellent agreement with those from the all-electron cases. For the platinum atom, the agreement between methods for the relative property differences between isomers are not as good as for the other atoms but are still quite reasonable. Tetrairidium Dodecacarbonyl. As a final example, we consider two isomers of a multinuclear third row transition metal complex, tetrairidium dodecacarbonyl.67 Figure 5a shows the QTAIM molecular graph of the Td symmetric form of this complex, in which the four iridium atoms are at the vertices of a tetrahedron and each is bonded to three carbonyl groups. Figure 5b shows the C3v isomer of this complex in which three of the carbonyl groups have migrated to bridging positions between Iridium atom pairs. The molecular graphs were calculated from electron densities from B3LYP/cc-pVDZ[C,O],53,54 cc-pVDZ-PP[Ir]68 equilibrium geometry, ECP-based wave functions augmented with 60-electron atomic core densities for each of the Ir atoms. Consistent with experimental observations of the Td tetrairidium dodecacarbonyl isomer,67 both the all-electron relativistic calculations and ECP calculations indicate (see Table 6) that the Td isomer is slightly more stable than the C3v isomer. Table 6 compares QTAIM atomic and bonding properties for the Td isomer using three different methods: all-electron,
ARTICLE
nonrelativistic B3LYP/cc-pVDZ[C,O], SARC-DKH[Ir];63 all-electron, relativistic B3LYP/cc-pVDZ[C,O], SARC-DKH[Ir]; and nonrelativistic B3LYP/cc-pVDZ[C,O],cc-pVDZPP[Ir] in which a 60-electron ECP was used for each Ir atom. In the last case, results from two methods of representing the Ir core electron densities are shown, namely Ir atomic core densities and tight core densities, eq 9, and they are essentially the same. As with the platin isomers, this ensures that the difference in QTAIM results between the ECP-based wave function and the all-electron wave functions is due only to the wave functions. Table 7 is similar to Table 6 except the data is for the C3v isomer of tetrairidium dodecacarbonyl while Table 8 shows the QTAIM atomic and bonding properties of the C3v bridging isomer of tetrairidium dodecacarbonyl relative to those in the Td isomer. Considering the slightly different valence basis sets on the Ir atoms in the all-electron vs ECP calculations as well as the difference between the nonrelativistic results and relativistic allelectron results, the agreement between the properties calculated using ECP-based wave functions augmented with Ir atomic core densities and the all-electron results is quite reasonable. As found for the cisplatin and transplatin examples, the absolute atomic kinetic energies are not comparable between methods due in part to basis set differences, but the relative atomic kinetic energies are comparable for the non-Ir atoms and are in reasonable agreement. Importantly, the relative QTAIM properties between the two isomers using ECP-based wave functions are consistently closer to those for the relativistic, all-electron wave functions than for the nonrelativistic ones.
’ CONCLUSIONS Electron densities of molecular wave functions from ab initio calculations using effective core potentials for heavy atoms can be problematic for QTAIM analyses due to the absence of electron density contributions from the heavy atom core electrons. The associated absence of strong nuclear maxima in the electron density distribution at the heavy atom nuclei results in spurious, low electron density critical points in the atomic core and sometimes spurious critical points even outside the atomic core. The missing electron density contribution from a given heavy atom’s ECP-modeled core electrons can be restored using an atomic core density obtained from an all-electron calculation for the free atom. The number of electrons comprising an atomic core can vary between ECP models, especially for very heavy elements. For an atom for which a small-core or medium-core ECP is used, simply representing the corresponding core electron density using a simple, tightly localized electron density function is actually sufficient to obtain a correct electron density topology and perform QTAIM analyses to obtain at least semiquantitatively meaningful results, but this is often not true when large-core ECPs are used. Relativistic B3LYP total atomic electron densities have been calculated for all elements of the periodic table and every subshell contribution to these total atomic densities has been separately fitted to a linear combination of spherical s-type Gaussian functions common to all subshells of a given element. This enables total atomic core densities to be systematically, consistently, and efficiently constructed and stored for any possible total atomic core defined in terms of closed subshells and 12892
dx.doi.org/10.1021/jp2040086 |J. Phys. Chem. A 2011, 115, 12879–12894
The Journal of Physical Chemistry A therefore be used to correct the electron densities obtained from ECP-based wave functions, for the vast majority of ECPs. The augmentation of ECP-based wave functions with appropriate atomic core densities has been automated in Gaussian24 for all elements of the periodic table and the use of these wave functions has been automated in the QTAIM analysis software package AIMAll.40 Results from QTAIM analyses of several heavy atom ECP-based molecular wave functions augmented with atomic core densities have been compared with corresponding results obtained by augmenting with simpler, tight core densities and also with QTAIM results calculated from allelectron nonrelativistic and relativistic wave functions. For small-core and medium-core ECPs, augmentation with atomic core densities gives QTAIM results in very close agreement with the QTAIM results obtained when augmenting with simplistic, tight core densities and generally gives good agreement with QTAIM results for all-electron wave functions using comparable valence basis sets, usually being within the range of differences between the all-electron nonrelativistic and relativistic results. For very heavy atoms, QTAIM results from ECP-based wave functions are generally closer to all-electron, relativistic QTAIM results than corresponding nonrelativistic ones. For large-core ECPs, agreement between QTAIM results from augmented ECP-based wave function and those from all-electron wave functions is less good but is usually at least semiquantitatively close. Comparison of QTAIM results from augmenting largecore ECP-based wave functions with atomic core densities versus tight core densities shows that the former generally gives better agreement with the all-electron results.
’ ASSOCIATED CONTENT
bS
Supporting Information. Table S1 lists spin multiplicities used for B3LYP/UGBS DKH2 calculations of free atoms of all elements of the periodic table. The corresponding KohnSham atomic orbitals were used for the atomic subshell electron density fittings. Document S2 describes the different schemes used for ensuring positive definiteness in the atomic subshell electron density fitting. Table S3 lists QTAIM atomic and bonding properties of transplatin using different methods. Table S4 lists QTAIM atomic and bonding properties of cisplatin using different methods. Table S5 lists QTAIM atomic and bonding properties of tetrairidium dodecacarbonyl (Td) using different methods. Table S6 lists QTAIM atomic and bonding properties of tetrairidium dodecacarbonyl (C3v) using different methods. Table S7 lists differences between QTAIM atomic and bonding properties of tetrairidium dodecacarbonyl (C3v) and tetrairidium dodecacarbonyl (Td) using different methods. This material is available free of charge via the Internet at http://pubs.acs.org
’ AUTHOR INFORMATION Corresponding Author
*E-mail:
[email protected].
’ REFERENCES (1) Pyykko, P. Chem. Rev. 1988, 88, 563. (2) Reiher, M.; Wolf, A. Relativistic Quantum Chemistry The Fundamental Theory of Molecular Science; Wiley VCH Verlag GmbH & Co. KGaA: Weinheim, 2009. (3) Dolg, M. In Effective Core Potentials. Modern Methods and Algorithms in Quantum Chemistry; Proceedings, 2nd ed.; Grotendorst, J.,
ARTICLE
Ed.; NIC Series, Proceedings, Vol. 3; John von Neumann Institute for Computing: Julich, 2000; ISBN 3-00-005834-6, p 507. (4) Cao, X.; Dolg, M. Wiley Interdiscip. Rev.: Comput. Mol. Sci. 2011, 1 (2), 200. (5) Bader, R. F. W. Atoms in Molecules A Quantum Theory; Oxford University Press: Oxford, U.K., 1990. (6) The Quantum Theory of Atoms in Molecules: From Solid State to DNA and Drug Design; Matta, C. F.; Boyd, R. J., Eds.; Wiley-VCH: Weinheim, 2007. (7) Srebrenik, S.; Bader, R. F. W. J. Chem. Phys. 1975, 63, 3945. (8) Bader, R. F. W.; Stephens, M. E. J. Am. Chem. Soc. 1975, 97, 7391. (9) Srebrenik, S.; Bader, R. F. W.; Nguyen-Dang, T. T. J. Chem. Phys. 1978, 68, 3667. (10) Bader, R. F. W.; Essen, H. J. Chem. Phys. 1984, 80, 1943. (11) Bader, R. F. W.; Nguyen-Dang, T. T. Adv. Quantum Chem. 1981, 14, 63. (12) Bader, R. F. W.; Nguyen-Dang, T. T.; Tal, Y. Rep. Prog. Phys. 1981, 44, 893. (13) Bader, R. F. W. Chem. Rev. 1991, 91, 893. (14) Fradera, X.; Austen, M. A.; Bader, R. F. W. J. Phys. Chem. A 1999, 103, 304. (15) Bader, R. F. W. Monatsch. Chem. 2005, 136, 819. (16) Cortes-Guzman, F.; Bader, R. F. W. Coord. Chem. Rev. 2005, 249, 633. (17) Ho, M.; Schmider, H.; Edgecombe, K. E.; Smith, V. H., Jr. Int. J. Quantum Chem. 1994, S28, 215. (18) Cioslowski, J.; Piskorz, P. Chem. Phys. Lett. 1996, 255, 315. (19) Cioslowski, J.; Piskorz, P.; Rez, P. J. Chem. Phys. 1997, 106 (9), 3607. (20) Bader, R. F. W.; Gillespie, R.; Martin, F. Chem. Phys. Lett. 1998, 290, 488. (21) Tiana, D.; Francisco, E.; Blanco, M. A.; Pendas, A. M. J. Phys. Chem. A 2009, 113, 7963. (22) Louie, S. G.; Froyen, S.; Cohen, M. L. Phys. Rev. B 1982, 26, 1738. (23) Delley, B. Phys. Rev. B 2002, 66, 155125. (24) 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.; Keith, T.; 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.; Parandekar, P. V.; Mayhall, N. J.; Daniels, A. D.; Farkas, O.; Foresman, J. B.; Ortiz, J. V.; Cioslowski, J.; , Fox, D. J. Gaussian Development Version; Gaussian, Inc.: Wallingford, CT, 2010. (25) Becke, A. D. J. Chem. Phys. 1993, 98, 5648–5652. (26) Lee, C.; Yang, W.; Parr, R. G. Phys. Rev. B 1988, 37, 785. (27) Vosko, S. H.; Wilk, L.; Nusair, M. Can. J. Phys. 1980, 58, 1200. (28) Stephens, P. J.; Devlin, F. J.; Chabalowski, C. F.; Frisch, M. J. J. Phys. Chem. 1994, 98, 11623. (29) de Castro, E. V. R.; Jorge, F. E. J. Chem. Phys. 1998, 108, 5225–29. (30) Douglas, M.; Kroll, N. M. Ann. Phys. (N Y) 1974, 82, 89–155. (31) Hess, B. A. Phys. Rev. A 1985, 32, 756–63. (32) Hess, B. A. Phys. Rev. A 1986, 33, 3742–48. (33) Jansen, G.; Hess, B. A. Phys. Rev. A 1989, 39, 6016. (34) Yang, R.; Rendell, A.; Frisch, M. J. J. Chem. Phys. 2007, 127, 074102. (35) Dunlap, B. I.; Connolly, J. W. D.; Sabin, J. R. J. Chem. Phys. 1979, 71, 4993. 12893
dx.doi.org/10.1021/jp2040086 |J. Phys. Chem. A 2011, 115, 12879–12894
The Journal of Physical Chemistry A
ARTICLE
(36) Dunlap, B. I. J. Chem. Phys. 1983, 78, 3140. (37) Dunlap, B. I. J. Mol. Struct. (THEOCHEM) 2000, 529, 37. (38) Dennington, Roy D.; Keith, Todd A.; Millam, John M. GaussView 5.0.9; Semichem, Inc.: Shawnee, KS, 2009. (39) Extended Wavefunction File Format (.wfx files): http://aim. tkgristmill.com/wfxformat.html. (40) AIMAll, Version 11.04.03; Keith, Todd A.; T. K. Gristmill Software: Overland Park, KS, 2011; aim.tkgristmill.com. (41) Bader, R; Popelier, P; Chang, C J. Mol. Struct. 1992, 255, 145. (42) Popelier, P. In Molecular Similarity in Drug Design; Dean, P. M., Ed.; Blackie: Glasgow, 1995; Vol. 215 (43) Popelier, P. L. A. J. Phys. Chem. A 1999, 103, 2883. (44) O’Brien, S. E.; Popelier, P. L. A. J. Chem. Inf. Comput. Sci. 2001, 41, 764. (45) O’Brien, S. E.; Popelier, P. L. A. Quantum Molecular Similarity: Use of Atoms in Molecules derived quantities as QSAR variables. European Congress on Computational Methods in applied Sciences and Engineering (ECCOMAS), Sep 2000, Barcelona, Spain (46) McLean, D.; Chandler, G. S. J. Chem. Phys. 1980, 72, 5639. (47) Raghavachari, K.; Binkley, J. S.; Seeger, R.; Pople, J. A. J. Chem. Phys. 1980, 72, 650. (48) Binning, R. C., Jr.; Curtiss, L. A. J. Comput. Chem. 1990, 11, 1206. (49) McGrath, M. P.; Radom, L. J. Chem. Phys. 1991, 94, 511–16. (50) Curtiss, L. A.; McGrath, M. P.; Blaudeau, J.-P.; Davis, N. E.; Binning, R. C., Jr.; Radom, L. J. Chem. Phys. 1995, 103, 6104. (51) Peterson, K. A. J. Chem. Phys. 2003, 119, 11099. (52) Martin, J. M. L.; Sundermann, A. J. Chem. Phys. 2001, 114, 3408. (53) Dunning, T. H., Jr. J. Chem. Phys. 1989, 90, 1007. (54) Balabanov, N. B.; Peterson, K. A. J. Chem. Phys. 2005, 123, 064107. (55) Hay, P. J.; Wadt, W. R. J. Chem. Phys. 1985, 82, 270. (56) Hay, P. J.; Wadt, W. R. J. Chem. Phys. 1985, 82, 299. (57) Stewart, C. A.; Harlow, R. L.; Arduengo, A. J. J. Am. Chem. Soc. 1985, 107, 5543. (58) Godbout, N.; Salahub, D. R.; Andzelm, J.; Wimmer, E. Can. J. Chem. 1992, 70, 560. (59) Sosa, C.; Andzelm, J.; Elkin, B. C.; Wimmer, E.; Dobbs, K. D.; Dixon, D. A. J. Phys. Chem. 1992, 96, 6630. (60) Peterson, K. A.; Figgen, D.; Goll, E.; Stoll, H.; Dolg, M. J. Chem. Phys. 2003, 119, 11113. (61) Alderden, R. A.; Hall, M. D.; Hambley, T. W. J. Chem. Educ. 2006, 83, 728. (62) Frisch, M. J.; Pople, J. A.; Binkley, J. S. J. Chem. Phys,. 1984, 80, 3265. (63) Pantazis, D. A.; Chen, X.-Y.; Landis, C. R.; Neese, F. J. Chem. Theory Comput 2008, 4, 908. (64) Roy, L. E.; Hay, P. J.; Martin, R. L. J. Chem. Theory Comput. 2008, 4, 1029. (65) Ehlers, A. W.; Bohme, M.; Dapprich, S.; Gobbi, A.; Hollwarth, A.; Jonas, V.; Kohler, K. F.; Stegmann, R.; Veldkamp, A.; Frenking, G. Chem. Phys. Lett. 1993, 208, 111. (66) Wadt, W. R.; Hay, P. J. J. Chem. Phys. 1985, 82, 284. (67) Holland, G. F.; Ellis, D. E.; Tyler, D. R.; Gray, H. B.; Trogler, W. C. J. Am. Chem. Soc. 1987, 109, 4276. (68) Figgen, D.; Peterson, K. A.; Dolg, M.; Stoll, H. J. Chem. Phys. 2009, 130, 164108.
12894
dx.doi.org/10.1021/jp2040086 |J. Phys. Chem. A 2011, 115, 12879–12894