Journal of Molecular Structure (Theochem) 535 (2001) 279±286
www.elsevier.nl/locate/theochem
The response of atomic electron densities to point perturbations in the external potential W. Langenaeker a,*, S. Liu b a
Eenheid Algemene Chemie, Vrije Universiteit Brussel, Faculteit Wetenschappen, Pleinlaan 2, 1050 Brussel, Belgium b Department of Chemistry, University of North Carolina, Chapel Hill, NC 27599-3290, USA Received 8 February 2000; revised 10 May 2000; accepted 10 May 2000
Abstract The response of the electron density to a perturbation in the external potential is studied for atoms. After the independence of calculational method is established, characteristics of the response function are described and connections with the polarizability and hardness are noted. q 2001 Elsevier Science B.V. All rights reserved. Keywords: Response function; Perturbation; Softness; Polarisability
1. Introduction Molecular chemical reactivity is often discussed in language implicating the electron density concept (for example, nucleophilicity, electrophilicity, electron donating and accepting character). The electronic charge distribution in a molecule therefore can be expected to explain, at least to some extent, its behavior in interaction with a reaction partner. This suggests that Density Functional Theory (DFT) [1], in which the electron density is the basic quantity, is particularly appropriate for discussing reactivity. A series of quantities, which are conveniently used when considering chemical reactivity, appear in a most natural way in DFT. First of all, a new theoretical basis is provided for the use of the frontier molecular orbitals (FMO) as reactivity indices. The frontier orbital concept was * Corresponding author. Janssen Research Foundation, Janssen Pharmaceutica NV, Turnhoutseweg 30, B-2340 Beerse, Belgium. Fax: 132-14-6025-57. E-mail address:
[email protected] (W. Langenaeker).
introduced by Fukui in the FMO theory [2,3]. The FMO are approximations to the Fukui functions which appear in DFT (vide infra). Secondly, a quantity appearing in density functional theory can be identi®ed with the hardness [4] of a system. The ideas of hardness and softness were introduced by Pearson in discussing acid±base reactions of the type A1 : B ! A : B where A is a Lewis acid and B is a Lewis base. Pearson proposed the HSAB (Hard and Soft Acids and Bases) principle, stating that hard acids prefer to react with hard bases, and soft acids prefer to react with soft bases. The knowledge of the hardness and softness of a system can help to explain and predict the chemical behavior of Lewis acids and bases. As will be shown below, knowledge of the linear response function v
r; r 0 ; de®ned as: dr
r 0 v
r; r
1 dn
r 0 N
0166-1280/01/$ - see front matter q 2001 Elsevier Science B.V. All rights reserved. PII: S 0166-128 0(00)00579-0
280
W. Langenaeker, S. Liu / Journal of Molecular Structure (Theochem) 535 (2001) 279±286
can be considered fundamental for the calculation of a number of reactivity descriptors. This function describes the response in the electron density r
r to a perturbation in the external potential n
r 0 ; and is essentially the screening function used in solid state physics [5]. The quantity v
r; r 0 is sometimes called the polarizability kernel. It is symmetric, v
r; r 0 v
r 0 ; r: As we are especially interested in reactivity, and as the perturbation in n
r associated with the presence of a reaction partner is not uniform, our approach to the evaluation of v
r; r 0 will be unconventional. A single point charge will be used to induce some change in n
r (vide infra) and a related change in r
r; thus mimicking a reaction partner. We hope to thereby obtain an insight in some local properties of v
r; r 0 for atoms via the calculation of Dr
r=Dn
r 0 N and to establish regularities in its characteristics. The dependence on the use of different calculational methods and basis sets will also be investigated. Before proceeding to the calculation of Dr
r=Dn
r 0 N ; we brie¯y summarize how reactivity indices arise in the frontier molecular orbital theory. In the simplest electron-repulsion-free MO theory, the total energy is expressed in terms of the Coulomb integrals Ar and the resonance integrals Brs [2]: XX X qr A r 1 2 prs Brs
2 E r
r,s
The in¯uence of a change in the Coulomb integral of atom r is given, up to second order, by: 2E 2E
3 DAr 1 DAr DAs DE 2Ar 2Ar 2As DE qr DAr 1
1 2
prs DAr DAs
4
where qr is the charge density on atom r and prs is the mutual atom±atom polarizability of atom r and s, which satis®es:
prs
2qr 2As
5
This formula describes the response in the charge density on atom r to a change in the coulomb integral of atom s. Thus, prs is essentially a special example of v
r; r 0 of Eq. (1). In the consideration of chemical reactivity in the
polarized state approximation (distinct, interacting systems), the in¯uence of a reagent R can be explained in terms of what one often calls a change in electronegativity of an atom considered. If an electrophilic reagent R 1 interacts with a substrate, its positive charge pulls the electrons of the system towards the atom it approaches (as if the electronegativity of the atom increases). This is described by the ®rst term in Eqs. (3) and (4). A second important factor is charge capacity, related to the inverse of the hardness [6], of the atom considered. The more easily the charge can be brought to the atom, the larger the change in q due to the given perturbation in A. The effect is the largest for the atom with the highest self polarizability prr : In DFT an analogue to Eq. (4), describing the change in E up to second order in Dn
r; is: Z dE Dn
r dr DEN; n
r dn
r N Z dE Dn
rDn
r 0 dr dr 0
6 1 dn
rdn
r 0 N DEN; n
r 1
Z
r
rDn
r dr
Z dr
r Dn
rDn
r 0 dr dr 0 dn
r 0 N
7
Note that the second term in Eq. (7) contains the linear response function (polarizability kernel) v
r; r 0 of Eq. (1). 2. Theory As mentioned in the foregoing, recently several new reactivity indices have been introduced [7,8], which ®nd their signi®cance within the HSAB-principle. Here we will establish the role of the response function v
r; r 0 as a fundamental quantity using two of these new quantities. First there is a local hardness related quantity, h(r) de®ned [7] in the isomorphic ensemble [9] as 1 dm
8 h
r N ds
r N with N the number of electrons, m the chemical potential and s
r the shape factor [10,11]. Second, there is
W. Langenaeker, S. Liu / Journal of Molecular Structure (Theochem) 535 (2001) 279±286 0
the quantity q
r; r de®ned as 1 dn
r q
r; r 0 N ds
r 0 N
9
This symmetric kernel is the inverse of the linear response function v
r; r 0 and will therefore be denoted as v21
r 0 ; r in the present paper. An unambiguous de®nition of local hardness, h
r; was provided by De Proft et al. [8] 2u
r h
r
10 2N s This leads to the formulas: Z h
r h
r 1 s
r 0
h
r 0 2 v 21
r; r 0 dr 0
11
and
h
r; r 0 h
r 2 h
r 0 2 v
r 0 ; r
12
where h
r; r 0 is the hardness kernel as earlier de®ned by [12]:
h
r; r 0
d2 F dr
rdr
r 0
13
with F the universal Hohenberg±Kohn functional. The softness-related counterparts of Eqs. (3) and (4) are given by [12]: Z
14 s
r s
r; r 0 dr 0 and s
r; r 0 v
r; r 0 2
s
rs
r 0 S
15
The conclusion is, based on Eqs. (8)±(15) that knowledge of v
r; r 0 is of fundamental importance for the effective calculation of all reactivity indices which ®nd their meaning in association with the HSAB principle. Calculations will enable us to get numerical veri®cation of the various relations, and also to get better insight in different aspects of reactivity. Furthermore, calculations will enable us to demonstrate the relation between the polarizability, a , of a system and the linear response function, a relation, which was recently established analytically within a DFT-framework [13]. In order to proceed with such calculations it is essential to better understand the quantity v
r; r 0 :
281
As a ®rst step, the change in r (r), Dr (r), upon the introduction of a small positive charge (0.01) mimicking some change in the external potential is investigated Dr
r rq
r 2 r0
r
16
where r0
r is the unperturbed density and rq
r is the perturbed density. The charge was taken to be as small as possible, in order to approximate the derivative as closely as possible, still generating meaningful numerical results. For atoms two cases will be considered: ®rst new charge is introduced at the nucleus of the atom, yielding information on Dr
r which is essentially due to intrinsic characteristics of atoms such as polarizability and hardness. The relation between Dr
r and the hardness also stresses the importance of the quantity v
r; r 0 for atoms as it is a measure for the maximum possible amount of the atom-to-atom electron drift in a molecule for a given atom. Secondly the new charge is placed at different distances from the nucleus of an atom, a study which is more directly relevant for considering reactivity, as the new charge can be regarded as mimicking the presence of a reaction partner. 3. Results and discussion All calculations were performed using DFT methods as implemented in the Gaussian 94 program [14]. The B3PW91 [15] functional was used in combination with the 6-311 1 G p [16] basis set throughout. The quality of the densities obtained by this method was ®rst established by means of a comparison of the obtained electronic structures with those of an earlier study using a series of basis sets and different postHartree±Fock methods [17]. For Li and B the B3PW91/6-311 1 G p Dr
r function was explicitly compared with results from other methods using different basis sets. First the quality of the B3PW91/6-311 1 G p densities is compared to that of the QCISD/6-311 1 G p [18] densities by means of a comparison of the position obtained for the core-valence border given [17], following Politzer and Parr [19], by the outermost minimum in the radial distribution function. This is the position of the single minimum for ®rst row atoms and the position of the second minimum for second
282
W. Langenaeker, S. Liu / Journal of Molecular Structure (Theochem) 535 (2001) 279±286
Table 1 Position of the core-valence border (a.u.)
Li Be B C N O F Ne Na Mg Al Si P S Cl Ar a b
QCISD
B3PW91
1.779 1.095 0.796 0.619 0.501 0.418 0.355
1.78 1.10 0.80 0.62 0.50 0.42 0.36 0.30
a
b
2.111 1.627 1.342 1.143 0.998 0.886 a
b
2.18 1.64 1.35 1.16 1.00 0.89 0.80
Not given. No second minimum observed.
row atoms. As is seen from Table 1, taking the difference in resolution of both results into account, essentially the same values are obtained with the different methods.
Next, the quality of the point-to-point information of Dr
r itself, calculated using different basis sets and methods, is investigated for the case of a pointcharge on the nucleus. It is seen from Fig. 1 that, as was observed for r
r [17], the in¯uence of the calculational method, whether Hartree±Fock (HF) [20,21], Moller±Plesset second order perturbation theory (MP2) [22] or B3PW91, is relatively small. The results of the B3PW91, HF and MP2 methods are very much the same. MP2-densities are known to be of a fair quality, so one concludes that the B3PW91densities, and especially the more method-sensitive differences in densities, are reliable when an appropriate basis set is used. In the past it has shown that this is so for the 6-311 1 G p basis set when looking at total densities. We seek to con®rm this for our density differences by looking at the basis set dependence using the B3PW91 functional. In Fig. 2 radial distribution function for Dr
r
4pr 2 Dr
r obtained using the B3PW91/6-31 1 G [16] and B3PW91/6-311 1 G p calculations are plotted for Li, and in Fig. 3 the B3PW91/aug-cc-pV5Z [23,24] and B3PW91/6311 1 G p results are plotted for B. There is an in¯uence when one uses different basis sets. Two remarks should be made, however. The
Fig. 1. Radial distribution function of Dr
r function for the Na atom calculated using HF, MP2 and B3PW91 methods in combination with 6311 1 G p basis set (distance between tick marks is 1 a.u.).
W. Langenaeker, S. Liu / Journal of Molecular Structure (Theochem) 535 (2001) 279±286
283
Fig. 2. Radial distribution function of Dr
r function for the Li atom calculated using B3PW91 method in combination with 6-311 1 G p and 631G basis set (distance between tick marks is 1 a.u.).
difference between the B3PW91/6-31G and the B3PW91/6-311 1 G p results is much larger than the difference between the B3PW91/aug-cc-pV5Z result, involving the use of an extremely large basis set, and the B3PW91/6-311 1 G p result. Furthermore, the effects of the basis set are concentrated mainly in
the core region of atoms, which is of incidental interest to us (vide infra), and mainly concerns differences in the absolute values of the extrema of the function and not the positions of extrema and other characteristic points of this function. We are especially interested in the region most
Fig. 3. Radial distribution function of Dr
r function for the B atom calculated using B3PW91 method in combination with 6-311 1 G p and aug-cc-pV5Z basis set (distance between tick marks is 1 a.u.).
284
W. Langenaeker, S. Liu / Journal of Molecular Structure (Theochem) 535 (2001) 279±286
Table 2 Hardness (eV), Polarisability (a.u.), and values (a.u.) for characteristical points in the function Dr
r (see Fig. 3) Atom
Hardness
Polarizability
a
x
y
Li Be B C N O F Ne Na Mg Al Si P S Cl Ar
3.01 4.9 4.29 6.27 7.3 7.54 10.41
149.05 42.767 17.571 8.888 5.059 2.800 1.667 1.054 163.418 73.729 39.614 22.906 14.530 9.966 6.976 5.14
2.04 1.3 0.94 0.72 0.58 0.48 0.42 0.36 2.74 1.97 1.66 1.44 1.26 1.12 1.02 0.92
1.53 0.95 0.67 0.49 0.39 0.31 0.25
4.13 2.57 1.99 1.63 1.37 1.19 1.03 0.93 4.29 3.19 2.89 2.53 2.19 1.95 1.73 1.59
a
a
2.85 3.75 3.23 4.77 5.62 6.22 8.3 a
2.17 1.57 1.29 1.11 0.97 0.85 0.73 0.69
Not given.
important for reactivity related/determining effects and so the characteristics of Dr
r we will focus on are mainly in the valence regions of atoms. The general behavior in the valence region of Dr
r for the case of a point-charge on the nucleus is the same for all atoms: an increase in electron density followed by a decrease when moving away from the nucleus. Values of the distance to the nucleus of the points x, y, and a (indicated in Fig. 3) are given in Table 2.
Comparison with the position of the core-valence border (Table 1) shows that in all cases this border is situated between the position of x and a, establishing the increase and decrease considered as a valence shell effect. One might expect a correlation between the characteristics of Dr
r; which is proportional to v
r; r 0 for a given perturbation Dn
r 0 ; and the polarizability, a , of the atom considered. The a -values for the atoms considered were calculated at the B3PW91/6311 1 G p level and are given in Table 2. Experimental values [1] for the global hardness are also included. One immediately notices the well-known correlation between the global hardness and the polarizability. When looking at the positions x, y, and a, one sees that the distance of these points to the nucleus decreases with decreasing polarisability. This would be expected on the basis of the behavior of the polarizability and the position of the core-valence border, each of which correlates well with the distances monitored for x, y and a, when moving through the periodic table. Note, however, that these geometrical parameters of Dr
r show a very good correlation with the a -values (Fig. 4). In general the effects of the perturbation in the external potential at the nucleus are closer to the nucleus as the polarizability decreases. This agrees with the idea of a smaller system being harder and thus less polarizable. If we look at the effect on r
r of a perturbation in the external potential at a certain distance from the
Fig. 4. Polarisability (a.u.) and value of distance ( p50) between position y and nucleus (a.u.) as a function of the atomic number.
W. Langenaeker, S. Liu / Journal of Molecular Structure (Theochem) 535 (2001) 279±286
285
Fig. 5. The Dr
r function for a point charge at different distances from the nucleus for the Li atom calculated using B3PW91 method in combination with 6-311 1 G p basis set (distance between tick marks is 1 a.u.).
nucleus, one can distinguish between two types of behavior. One essentially shows the characteristics of the charge on the nucleus case (boldface line in Fig. 5), typical for the charge placed in the core region of the atom. The other essentially shows a increase in density on the side of the charge, compensated by a decrease on the opposite side of the nucleus (thin lines in Fig. 5), typical for a charge placed in the valence region of an atom. It is interesting that the shape of the Dr
r-function is insensitive to the distance at which the charge is placed. The positions of maxima and minima stay essentially the same. (Charges placed at a distance of up to 6 Bohr were considered in this study. These results were not included in Fig. 5, as the large distance in combination with the small charge generates very small values.) The same behavior was observed earlier in a study of the Fermi hole, where the shape was found to be insensitive to changes within one shell in the position of the probe electron [25]. Another remarkable resemblance between these two cases is the absence of a cusp at the position of the added charge or the probe electron. Furthermore, the position of the maximum on the side of the charge coincides with the position of the maximum detected at point a in the case of a charge at the nucleus. Another characteristic of Dr
r is the integrated
charge delocalization associated with a perturbation in the external potential. This was not to be checked in the present study, for two reasons. As was shown in the study of the quality of calculated Dr
r functions, the x, y and a positions are relatively insensitive to method and basis set effects. The absolute values of the function however is sensitive to basis set effects and so a more elaborate study of these effects is needed before conclusions are drawn. An additional reason for the case of the charge on the nucleus resides in the relatively large changes it generates in the core region of the atom. 4. Conclusions The point-to-point information about Dr
r; has been found to be relatively insensitive to the in¯uence of the calculational method. In the dependence on basis set, the general shape and positions of characteristic point (such as positions of extrema) in Dr
r are reproduced using the moderately sized 6-311 1 G p basis set. However, if one is interested in accurate absolute values of the extrema of Dr
r a more extensive basis set is needed. When one applies a perturbation in the external
286
W. Langenaeker, S. Liu / Journal of Molecular Structure (Theochem) 535 (2001) 279±286
potential at the nucleus of an atom, the general behavior in the valence region of the Dr
r functions is found to be the same for all atoms: an increase in electron density near the nucleus and a decrease of density at large distance from the nucleus. Values of the distance to the nucleus of some characteristic points in the valence region correlate very well with polarizability values. When looking at the effect on r
r of a perturbation in the external potential at a certain distance from the nucleus, two types of results are found. When the charge is placed in the core region of the atom, the behavior is like the charge on nucleus case. When the charge is placed in the valence region of the atom, the effects are relatively insensitive to precisely where the charge is placed. Acknowledgements W.L. wishes to acknowledge the Fund for Scienti®c Research-Flanders (Belgium) (F.W.O.) for a position as postdoctoral fellow and was additionally supported by Fulbright and NATO travel grants enabling a stay at the University of North Carolina at Chapel Hill. Support is acknowledged of a grant from the Petroleum Research Fund of the American Chemical Society and the North Carolina Supercomputing Center. The author wishes to thank Professor R.G. Parr for numerous helpful discussions. References [1] R.G. Parr, W. Yang, Density Functional Theory of Atoms and Molecules, Oxford University, Oxford, 1989.
[2] I. Fleming, Frontier Orbitals and Organic Reactions, Wiley, London, 1976. [3] K. Fukui, T. Yonezawa, H. Shingu, J. Chem. Phys. 20 (1952) 722. [4] R.G. Pearson, J. Am. Chem. Soc. 85 (1963) 3533. [5] N.W. Ashcroft, N.D. Mermin, Solid State Physics, Saunders College, Philadelphia, 1976. [6] P. Politzer, J.S. Murray, M.E. Grice, in: K.D. Sen (Ed.), Chemical Hardness, Springer, Berlin, 1993, p. 101. [7] B.G. Baekelandt, A. Cedillo, R.G. Parr, J. Chem. Phys. 103 (1995) 8548. [8] F. De Proft, S. Liu, R.G. Parr, J. Chem. Phys. 107 (1997) 3000. [9] R.F. Nalewajski, R.G. Parr, J. Chem. Phys. 77 (1982) 399. [10] R.G. Parr, L. Bartolotti, J. Phys. Chem. 87 (1983) 2810. [11] A. Cedillo, Int. J. Quant. Chem. Symp. 28 (1994) 231. [12] M. Berkowitz, R.G. Parr, J. Chem. Phys. 88 (1988) 2554. [13] W. Langenaeker, F. De Proft, F. Tielens, P. Geerlings, Chem. Phys. Lett. 288 (1998) 628. [14] M.J. Frisch, G.W. Trucks, H.B. Schlegel, P.M.W. Gill, B.G. Johnson, M.A. Robb, J.R. Cheeseman, T. Keith, G.A. Petersson, J.A. Montgomery, K. Raghavachari, M.A. Al-Laham, V.G. Zakrzewski, J.V. Ortiz, J.B. Foresman, J. Cioslowski, B.B. Stefanov, A. Nanayakkara, M. Challacombe, C.Y. Peng, P.Y. Ayala, W. Chen, M.W. Wong, J.L. Andres, E.S. Replogle, R. Gomperts, R.L. Martin, D.J. Fox, J.S. Binkley, D.J. Defrees, J. Baker, J.P. Stewart, M. Head-Gordon, C. Gonzalez, J.A. Pople, Gaussian 94, Revision C.3, Gaussian Inc., Pittsburgh PA, 1995. [15] A.D. Becke, J. Chem. Phys. 98 (1993) 5648. [16] W.J. Hehre, L. Radom, P.v.R. Schleyer, J.A. Pople, Ab Initio Molecular Orbital Theory, Wiley, New York, 1986. [17] F. De Proft, P. Geerlings, Chem. Phys. Lett. 220 (1994) 405. [18] J.A. Pople, M. Head-Gordon, K. Raghavachari, J. Chem. Phys. 87 (1987) 5968. [19] P. Politzer, R.G. Parr, J. Chem. Phys. 64 (1976) 4634. [20] C.C.J. Roothan, Rev. Mod. Phys. 23 (1951) 69. [21] J.A. Pople, R.K. Nesbet, J. Chem. Phys. 22 (1959) 571. [22] C. Moller, M.S. Plesset, Phys. Rev. 46 (1934) 618. [23] D.E. Woon, T. Dunning, J. Chem. Phys. 100 (1994) 2975. [24] R.A. Kendall, T.H. Dunning Jr., R.J. Harrison, J. Chem. Phys. 96 (1992) 6796. [25] W.L. Luken, D.N. Beratan, Theor. Chim. Acta 61 (1982) 265.