How far do electrons delocalize?

How far do electrons delocalize? Benjamin G. Janesko, Giovanni Scalmani, and Michael J. Frisch Citation: The Journal of Chemical Physics 141, 144104 (2014); doi: 10.1063/1.4897264 View online: http://dx.doi.org/10.1063/1.4897264 View Table of Contents: http://scitation.aip.org/content/aip/journal/jcp/141/14?ver=pdfcov Published by the AIP Publishing Articles you may be interested in How accurate is the strongly orthogonal geminal theory in predicting excitation energies? Comparison of the extended random phase approximation and the linear response theory approaches J. Chem. Phys. 140, 014101 (2014); 10.1063/1.4855275 An accurate first principles study of the geometric and electronic structure of B 2 , B 2 − , B 3 , B 3 − , and B 3 H : Ground and excited states J. Chem. Phys. 132, 164307 (2010); 10.1063/1.3389133 Improved version of a local contracted configuration interaction of singles and doubles with partial inclusion of triples and quadruples J. Chem. Phys. 132, 034108 (2010); 10.1063/1.3292605 Electronic structure of the trimethylenemethane diradical in its ground and electronically excited states:Bonding, equilibrium geometries, and vibrational frequencies J. Chem. Phys. 118, 6874 (2003); 10.1063/1.1561052 Low lying electronic states of rare gas–oxide anions: Photoelectron spectroscopy of complexes of O − with Ar, Kr, Xe, and N 2 J. Chem. Phys. 117, 2619 (2002); 10.1063/1.1491410

This article is copyrighted as indicated in the article. Reuse of AIP content is subject to the terms at: http://scitation.aip.org/termsconditions. Downloaded to IP: 193.0.65.67 On: Thu, 04 Dec 2014 10:56:27

THE JOURNAL OF CHEMICAL PHYSICS 141, 144104 (2014)

How far do electrons delocalize? Benjamin G. Janesko,1,a) Giovanni Scalmani,2 and Michael J. Frisch2 1 2

Department of Chemistry, Texas Christian University Fort Worth, Texas 76129, USA Gaussian, Inc., 340 Quinnipiac St., Bldg. 40 , Wallingford, Connecticut 06492, USA

(Received 11 June 2014; accepted 23 September 2014; published online 10 October 2014) Electron delocalization is central to chemical bonding, but it is also a fundamentally nonclassical and nonintuitive quantum mechanical phenomenon. Tools to quantify and visualize electron delocalization help to understand, teach, and predict chemical reactivity. We develop a new approach to quantify and visualize electron delocalization in real space. Our electron delocalization range function EDR(r ; u) quantifies the degree to which electrons at point r in a calculated wavefunction delocalize over length scale u. Its predictions are physically reasonable. For example, EDR(r ; u = 0.25 bohr) is close to one at points r in the cores of first-row atoms, consistent with the localization of core electrons to ∼0.25 bohr. EDR(r ; u = 1 bohr) is close to one at points r in typical covalent bonds, consistent with electrons delocalizing over the length of the bond. Our approach provides a rich representation of atomic shell structure; covalent and ionic bonding; the delocalization of excited states, defects, and solvated electrons; metallic and insulating systems; and bond stretching and strong correlation. © 2014 AIP Publishing LLC. [http://dx.doi.org/10.1063/1.4897264] I. INTRODUCTION

Electron delocalization is central to our understanding of chemistry. Covalent bonds form when electrons delocalize between atoms. Valence electrons are more delocalized than core electrons. However, our intuition about the world is based on macroscopic objects that do not delocalize. “Delocalization” is therefore a nonclassical and nonintuitive phenomenon. This fundamentally nonintuitive aspect of chemistry’s most fundamental concept complicates our efforts to teach chemical physics, to understand novel chemical bonding situations, and to improve electronic structure approximations. Several theoretical tools have been developed to build insight into electron delocalization.1, 2 These tools range from Lewis structures,3 to molecular orbital and natural bond orbital theory,4 to descriptors including the Fukui function,5, 6 the electron localization function (ELF),7–10 and the quantum theory of atoms in molecules (QTAIM).11 However, even with these tools, it is difficult to answer the question “How far do electrons delocalize?” We present a direct real-space measure of electron delocalization. Our electron delocalization range function EDR(r ; u) quantifies the degree to which an electron at point r delocalizes over length scale u. Plots and averages over EDR(r ; u) evaluated from calculated wavefunctions illustrate familiar concepts, and provide new insights into delocalization and chemical bonding. II. THEORY A. Background

(r1 , r2 . . . rN ).12 In particular, the nonlocal one-particle density matrix quantifies the probability that any one of a system’s indistinguishable electrons delocalizes between points r and r : γ (r , r ) ≡ N d 3 r2 . . . d 3 rN (r , r2 . . . rN ) ∗ (r , r2 . . . rN ). (1) Diagonal elements limr →r γ (r , r ) give the electron density ρ(r ), i.e., the classical probability density for finding an electron at r. Off-diagonal elements quantify electron delocalization. Points where r and r are on different atoms can define chemical bond orders: γ (r , r ) > 0 is a bonding interaction, and γ (r , r ) < 0 is an antibonding interaction.13 The exchange hole hX (r , r ) = −|γ (r , r )|2 /ρ(r ) provides a quantum-mechanical correction to the pair probability of finding one electron at r and a second electron at r .14, 15 Unlike molecular orbitals, the one-particle density matrix is uniquely defined in many-electron systems.4 However, the density matrix is still a nonclassical and nonintuitive function of six variables. B. Quantifying delocalization in real space

We define the electron delocalization range EDR(r ; u) to quantify the extent to which an electron at point r delocalizes over length scale u. We obtain it by contracting the density matrix γ (r , r ) with a test function of the distance of delocalization |r − r |: (2) EDR(r ; u) = d 3 r gu (r , r )γ (r , r),

Quantum mechanics postulates that all properties of an N-electron system can be obtained from its wavefunction

gu (r , r ) ≡

a) [email protected]

0021-9606/2014/141(14)/144104/13/$30.00

141, 144104-1

2 π u2

3/4 ρ

−1/2

|r − r |2 . (3) (r ) exp − u2 © 2014 AIP Publishing LLC


144104-2

Janesko, Scalmani, and Frisch

J. Chem. Phys. 141, 144104 (2014)

FIG. 1. EDR(r ; u) and ELF plotted at points r in the atomic plane of naphthalene. Light colors correspond to values near 1, dark colors to values near 0. The N-electron wavefunction is calculated with standard PBE0/6-311G(2d,2p) DFT. The roughness of the contours is due to the finite grid of points and the limits of the interpolation.

Length scale u sets the test function’s Gaussian decay in |r − r |. The prefactors in Eq. (3) ensure that |EDR(r ; u)|2 ≤ 1. Reference 16 and Secs. SI-II and SI-III of the supplementary material17 prove this inequality and provide technical details of our analytic integration over r . Briefly, if a system’s orbitals are expanded in standard atom-centered (AO) basis functions,18, 19 the integral in Eq. (2) simply requires the overlaps between the AOs and Gaussian functions centered at each gridpoint {ri }. We evaluate these integrals analytically using the standard PRISM algorithm.20 Plotting EDR(r ; u) at representative u in naphthalene (Figure 1) shows how it illustrates and quantifies electron delocalization and chemical bonding. EDR(r ; u = 0.25 bohr) (Figure 1(a)) is close to one in the carbon atom cores and close to zero everywhere else. This is consistent with the localization of carbon’s core electrons to short length scales u ∼ 0.25 bohr.21 EDR(r , u = 1.0 bohr) (Figure 1(b)) is close to zero in the carbon atom cores, and close to one in the C–C and C–H bonds. This is consistent with the delocalization of valence electrons over ∼1 bohr, approximately the

length of a typical chemical bond. Electrons on the periphery of the molecule (Figure 1(c)) delocalize over 2 bohrs or more. The peaks in EDR (r ; u = 2.0 bohrs) at the ring centers rcenter , where ρ(rcenter ) is small, appear to occur because γ (rcenter , r ) is dominated by points r near the C–C bonds. These bonds’ midpoints are located ∼2 bohrs from rcenter . This leads to a non-negligible contribution to the integrand of Eq. (2) at larger u, combined with a large contribution from the ρ −1/2 (rcenter ) prefactor of Eq. (3). Figure SI-1 of the supplementary material shows naphthalene’s calculated geometry.17 Figure SI-2 of the supplementary material illustrates the ELF and EDR along the axes of naphthalene’s shorter and longer C–C bonds.17 Both the ELF and EDR show small differences between the two bonds. Weighted averages over the EDR can quantify the average number of electrons that delocalize over length scale u: EDR(u) =

d 3 r ρ(r ) EDR(r ; u).

(4)


144104-3


These averages are particularly useful for comparing related systems. We define EDR(A, B; u) as the difference between the EDR(u) of two related systems A and B. If the electrons in system A are on average more localized than those in B, EDR(A, B; u) will have a positive peak at short length scales u and a negative peak at long u.

C. Relation to other descriptors

The EDR is a nontrivial extension of existing descriptors of electron delocalization. The vast literature on these descriptors has been recently reviewed.1, 2 For our purposes, these descriptors may be broadly grouped into two classes. The first class of descriptors are constructed from functions of one position r: the electron density ρ(r ), its gradients with respect to position or electron number, and the kinetic energy density τ (r ). Examples of this class of descriptors include QTAIM basins and attractors,11 the Fukui function,5, 6 the noncovalent interaction index,22 the density inhomogeneity,23 the localized orbital locator,24 and the ELF.7–10 Kohout’s electron localizability index,2, 25 which yields identical values to the ELF at the Hartree-Fock level of calculation, is constructed from the full pair density and thus is arguably not a member of this first class. The EDR may be viewed as generalizing this first class of descriptors. The “ingredients” of this class of descriptors, i.e., ρ(r ), τ (r ), and so on, can be constructed by replacing gu (r , r ) in Eq. (2) with nearly diagonal test functions. To illustrate, test function δ(r − r ) recovers the electron density ρ(r ). Test function (1/2)δ(r − r )∇r · ∇r recovers τ (r ).8 Our test function instead includes explicit information about delocalization along |r − r |, information unavailable to this first class of descriptors. The value of explicit information about nonlocality is illustrated by Figure 1(d), which shows the ELF in naphthalene. In contrast to the EDR, the ELF is close to 1 in all three regions of interest: C atom cores, C–C bonds, and C–H bonds. Knowing only that “the ELF is near 1 at point r does not tell whether point r is in an atomic core, a C–H bond, or some other region of localized electron density. In contrast, knowing that “EDR(r ; u = 1.0 bohr) is near 1 at point r indicates that point r is in a valence region, not in the core of a first-row atom. While one may use chemical intuition to distinguish the core and valence regions in Figure 1(d), we will suggest that the additional explicit information about delocalization provided by the EDR can help build new insights into more complicated systems. The second class of descriptors are constructed from functions of two positions {r , r }: the exchange hole hX (r , r ), the electron pair density, and related quantities. Examples of this class of descriptors include functions of three variables generated by choosing a reference position r (e.g., Ref. 15) or by averaging over points r ∈ A assigned11 to atom A (domain-averaged Fermi hole).26–28 The second class of descriptors also include scalars that quantify A-B bonds by averaging over r ∈ A, r ∈ B (localization and delocalization indices29–32 ). The EDR complements these descriptors by replacing the choice of special points or atoms with a choice

J. Chem. Phys. 141, 144104 (2014)

of delocalization length scale u. The EDR also has the benefit of being constructed entirely from one-electron quantities, aiding its evaluation in large systems. Schmider’s parity function33 is a particularly relevant member of this second class of descriptors. This function is equal to Eq. (2) with a constant test function g(r , r ) = 1 and a change of coordinate γ (r + (r − r)/2, r − (r − r)/2). We suggest that the EDR incorporates some of the physical insights of the parity function, while avoiding some of the challenges inherent to interpreting the parity function.34 Another relevant member of this second class of descriptors, proposed by Proynov and co-workers,35 is based on the effective exchange-hole normalization from Becke’s 2005 model of nondynamical correlation.36–38 This model is conexact exstructed from γ (r , r ) via the conventional-gauge change energy density eX (r ) = −(1/2) d 3 r |γ (r , r )|2 |r − r |−1 . It projects the exact exchange hole onto a semilocal model,39 quantifying delocalization past a the length scale set by the semilocal model.36, 37 This projection also inspired our “density matrix similarity metric,” which projects γ (r , r ) (rather than the exchange hole) onto a semilocal model density matrix.40, 41 The EDR complements these descriptors by scanning over multiple delocalization length scales u, rather than choosing a single length scale defined by a model exchange hole.

III. COMPUTATIONAL DETAILS

The remainder of this work illustrates the EDR in systems of fundamental interest for chemistry and physics. Calculations use the development version of the GAUSSIAN suite of programs.42 All calculations invoke the Born-Oppenheimer approximation. No corrections are included for relativistic or continuum solvent effects. One-electron orbitals are expanded in standard Pople,18 correlation consistent,19 and UGBS43 AO basis sets. Most calculations use single-determinant wavefunctions from Hartree-Fock theory or generalized KohnSham44–46 density functional theory (DFT). DFT calculations use the LSDA,47 PBE,48 PBE0,48–50 or ωB97XD51 approximate exchange-correlation functionals.52 Some calculations use nearly exact (full configuration interaction, FCI) wavefunctions. The one-particle density matrix from correlated and excited-state calculations is obtained from the Z-vector method.53 The ELF and EDR are evaluated at points {ri } on standard DFT numerical integration grids.54 The ELF is a function of three independent variables ˆ y, ˆ zˆ } are unit vectors in real r = x xˆ + y yˆ + z zˆ , where {x, space. EDR(r ; u) is a function of four variables: x, y, z, and the delocalization range u. 2D contour plots of these functions can vary at most two independent variables. We present two classes of contour plots based on two different choices of independent variables. The first class are exemplified by Figure 1 and by the ELF plots in Ref. 8. These plot the ELF, and the EDR at fixed u, at points {x, y} lying in a plane. The second class of contour plot is used exclusively for EDR(x; u). The abscissa is position x along a high-symmetry coordinate such as a chemical bond axis. The ordinate is length scale u. To aid


144104-4


J. Chem. Phys. 141, 144104 (2014)

FIG. 2. Kr atom. (a) Contour plot of EDR(r; u) plotted as a function of distance r from the nucleus and delocalization length scale u. Light colors correspond to EDR(r; u) near one, dark colors correspond to EDR(r; u) near zero. (b) ELF plotted as a function of r, and EDR(r; u) at the four u values shown as horizontal lines in (a). N-electron wavefunctions from HF/cc-pV5Z all-electron calculations.

readability, we plot the first class of contour plot in blue and the second class in green.

IV. RESULTS A. Atomic shell structure

Figure 2 illustrates the ELF and EDR in a spherically symmetric ground-state krypton atom. Figure 2(a) shows a two-dimensional contour plot of EDR(r; u), where r is the distance from the nucleus and u is the delocalization length scale. The plot shows four maxima corresponding to krypton’s four electronic shells. The maxima occur approximately along the diagonal of the plot, at increasing distances from the nucleus r = 0, 0.11, 0.44, 1.4 bohrs and at increasing delocalization lengths u = 0.04, 0.12, 0.34, 1.7 bohrs. (Figure 2(a) shows these u values as horizontal lines.) EDR(r; u) thus illustrates and quantifies the increasing delocalization of shells further from the nucleus. Figures SI-3– SI-8 of the supplementary material show similar results for He, Be, Ne, Ar, Zn, and Xe atoms.17 Section SI-VI of the supplementary material argues that the EDR(r; u) < 0 regions seen in some of these plots arise from filled atomic shells’ “antibonding” interaction.17 Figure 2(b) shows a one-dimensional plot of the ELF versus distance r from the nucleus. Like the EDR(r; u) in Figure 2(a), the ELF has four peaks corresponding to krypton’s four atomic shells.7 However, the EDR(r, u) in Figure 2(a) provides the additional quantitative insight that the first shell delocalizes over 0.04 bohr, the second shell delocalizes over 0.12 bohrs, and so on. Figure 2(b) also highlights the EDR’s ability to focus on electrons delocalized over a particular region of interest, by plotting the four “slices” through the maxima in EDR(r; u) shown in Figure 2(a). Each of these “slices” approximately corresponds to one of the four peaks in the ELF. Expectation values EDR(u) from the four maxima in EDR(r; u) (Figure 2(a)) are roughly consistent with the expected occupancies of krypton’s four shells. The valence n = 4 shell, which peaks at u = 1.7 bohrs, has EDR(u

= 1.7 bohrs) = 4.8. This is consistent with krypton’s valence shell containing four ↑-spin electrons. The n = 3 shell has EDR(u = 0.34 bohr) = 9.5, consistent with krypton’s n = 3 shell containing nine ↑-spin electrons. The agreement is degraded for the n = 1 and n = 2 shells: EDR(u = 0.12 bohr) = 6.9 and EDR(u = 0.04 bohr) = 2.9 are somewhat different from the expected occupancies of four and one ↑-spin electrons. This is consistent with Figure 2(b), which shows that EDR(r; u = 0.12 bohr) has significant contributions from the n = 1 and n = 3 shells, while EDR(r; u = 0.04 bohr has a significant contribution from the n = 2 shell. We are exploring other test functions in Eq. (3) to provide a better overlap with atomic cores. Figure 3 illustrates some of the chemical insight available from the EDR. The figure plots EDR(r; u) for zinc atom following the conventions of Figure 2(a). Comparing Figures 2(a) and 3 shows that the four maxima in zinc’s EDR are shifted rightward to larger distances from the nucleus, and upward to longer delocalization lengths. Zinc’s valence shell peaks at {r = 2.6 bohrs, u = 4.7 bohrs}, significantly more delocalized than krypton’s valence shell at {r = 1.4 bohrs,

FIG. 3. EDR(r; u) in Zn atom, details as in Figure 2(a).


144104-5


J. Chem. Phys. 141, 144104 (2014)

FIG. 4. EDR(r ; u) and ELF in the Kevan structure for a solvated electron (H2 O)− 6 . Results are plotted at points r in a plane containing four O–H bonds. The three u values correspond to maxima in EDR(x; u) along the O–H · · · H–O axis (Figure 5(b)).

u = 1.7 bohrs}. This illustrates that zinc’s valence electrons are more delocalized than krypton, consistent with its facile ionization to Zn2 + . B. Delocalized and solvated electrons

The EDR shows particular promise for visualizing regions of delocalized electron density. Such regions can occur in a variety of systems, including solvated electrons,55 electrides,56–58 and defect sites in solids (Figure 12). Figure 4 illustrates this ability for the Kevan59 model of a hydrated electron (H2 O)− 6 . This model treats an octahedral water cluster with six O–H bonds pointed towards the cluster center. Figure SI-10 of the supplementary material shows the calculated structure.17 Wavefunctions are evaluated using DFT with the LC-ωPBE functional.60, 61 Geometries and basis sets are taken from Ref. 57. Figures 4(a)–4(c) show the EDR for majority-spin electrons in a plane containing four O–H bonds. Results are plotted at three different delocalization ranges u=0.26, 1.2, 3.8 bohrs. These three u values correspond to the three sets of maxima in EDR(x, u) seen for points x along a high-symmetry O–H · · · H–O axis (Figure 5). EDR(r ; u = 0.26 bohr)

peaks in the oxygen atom cores (Figure 4(a)). EDR (r ; u = 1.2 bohrs) peaks in the valence regions (Figure 4(b)). EDR (r ; u = 3.8 bohrs) gives a single peak in the center of the structure, corresponding to the delocalized electron density (Figure 4(c)). As evidence that this is not an artifact of the approximation, Figure SI-11 of the supplementary material17 shows that the ELF and EDR are both negligible in the center of neutral (H2 O)6 , where the solvated electron is removed. In contrast to the EDR, the ELF (Figure 4(d)) shows core, valence, and delocalized regions simultaneously. The EDR complements the ELF by visualizing and quantifying each region’s degree of delocalization. To illustrate, Figure 5 shows the EDR in stretched and compressed Kevan structures, obtained by rigidly shifting the H2 O groups towards or away from the structure center. The figure plots EDR(x; u) on points x along an O–H · · · H–O axis intersecting the structure center. The O atom core electrons are visible at the bottom of the plots at small u. The valence regions are at the left and right of the plots at modest u. The solvated electron is visible in the structure center at large u. Compressing the Kevan structure by 0.5 Å (Figure 5(a)) localizes the solvated electron, moving the EDR maximum at the structure center from u = 3.8 bohrs (Figure 5(b)) to u = 2.7 bohrs. Expanding the structure


144104-6


J. Chem. Phys. 141, 144104 (2014)

FIG. 5. EDR(x; u) in compressed and expanded Kevan structures. Abscissa x is the position along a high-symmetry O–H · · · H–O axis, ordinate u is the delocalization range. The H2 O molecules are shifted from their equilibrium positions, relative to the structure center, by (a) −0.5, (b) 0, (c) +0.5, and (d) 1.0 Å. Horizontal lines denote the peak positions.

by 0.5 or 1.0 Å (Figures 5(c)–5(d)) delocalizes the solvated electron, moving the maximum to u = 4.3 bohrs and u = 6.3 bohrs, respectively. Figure SI-12 of the supplementary material shows the ELF along these O–H · · · H–O axes in stretched and compressed Kevan structures.17 The ELF shows a peak in the structure center, as expected. Figures SI-13 and SI-14 of the supplementary material show that the EDR and ELF are negligible in the centers of stretched and expanded neutral (H2 O)6 .17

Figure 6(c) plots the difference EDR (boron nitride, graphene; u). The plot shows a positive peak at small delocalization lengths u arising from the localized, insulating boron nitride electrons; and a negative peak at large u from the delocalized graphene electrons. This further illustrates the EDR’s ability to quantify electron delocalization. Figures SI-15– SI17 of the supplementary material illustrate calculated geometries, the ELF, and the EDR at other length scales u.17

C. Metals and insulators

D. Aromaticity and π bonding

Figure 6 illustrates how electrons in graphene are more delocalized than electrons in insulating boron nitride. Figures 6(a) and 6(b) compare EDR (r ; u = 1.1 bohrs) in a single unit cell of isolated, infinite, periodic sheets of graphene and boron nitride.62 EDR (r ; u = 1.1 bohrs) shows the C, N, and B atom cores as minima, and shows the bonding electrons as maxima along the C–C and B–N bond axes. These maxima occur at the C–C bond midpoints, and at points relatively close to the nitrogen atom in B–N bonds. This is consistent with the fact that the B–N bonds are more ionic than C–C bonds.

Much of the research on electron delocalization focuses on quantifying the concept of aromaticity.1, 63–65 Figure 7 illustrates the ELF and EDR in two model systems: aromatic D6h benzene, and a fictitious D3h “cyclohexatriene” with C–C and C=C bonds fixed to 1.54 and 1.36 Å.15, 63 Figure SI-18 of the supplementary material shows the calculated geometries.17 Figures 7(a)–7(c) show contour plots of EDR(x; u) along the axes of the benzene C–C bond, the cyclohexatriene C=C bond, and the cyclohexatriene C–C bond. In contrast to Figure 2(a), the horizontal lines are included solely to improve


144104-7


J. Chem. Phys. 141, 144104 (2014)

FIG. 6. EDR(r , u = 1.1 bohrs) for graphene (a) and boron nitride (b) sheets. Results are plotted for points r in the atomic plane, in a single unit cell of the infinite, isolated sheet. (c) EDR(boronnitride, graphene; u) plotted as a function of delocalization range u. Wavefunctions evaluated with PBE/6-31G DFT.

readability. Core electrons are visible at small u. Figures 7(d)– 7(f) show the ELF along the bond axes. Both the ELF and EDR suggest that the cyclohexatriene “C=C” bond is comparable to the benzene C–C bond. In contrast, the cyclohexatriene “C–C” bond is more delocalized, with the EDR peak around the bond midpoint shifted slightly upward to larger u. The presence of substantial delocalized electron density in the putative “C–C” bonds is consistent with previous reports of significant π electron delocalization in D3h cyclohexatriene.63 The ELF does not directly quantify this increased delocalization. Previous studies have shown that such systems require rather sophisticated topological analyses of the ELF.64 Figures SI-19 and SI-20 of the supplementary material17 show the ELF and the EDR at representative u, plotted in the nuclear plane and 0.5 bohr above the nuclear plane, where the benzene π density attains its maximum value.15 The two planes are qualitatively similar in the valence region, consistent with the nonuniqueness of the σ and π representation of double bonds.66 Figures SI-21 of the supplementary material17 illustrates EDR(benzene, cyclohexatriene; u), showing that benzene’s electronic structure is overall more localized than the fictitious and unstable distorted cyclohexatriene. Figure SI-22 of the supplementary material17 illustrates EDR(phenanthrene, anthracene; u), showing that the more stable, arguably more aromatic65 phenanthrene is on average more localized.

E. Bond dissociation and strong correlation

The EDR also illustrates and quantifies strong correlation in chemical bonds. Strong correlation is central to modern electronic structure theory.67–69 The H–H bond in stretched, singlet, D∞h -symmetric H2 provides a canonical illustration of strong correlation. References 27, 36, 67 and 70–72 exemplify the enormous literature on this model system. We begin by briefly reviewing previous results. At the equilibrium H–H bond length, the two electrons in H2 delocalize across both atoms. In the exact wavefunction, increasing the H–H distance breaks the bond, re-localizes the oneparticle density matrix, and gives an average of 1/2 ↑-spin and 1/2 ↓-spin electrons on each atom. Left-right correlation

between the two electrons shows up in the opposite-spin twoparticle density matrix. Simple restricted Hartree-Fock (RHF) theory calculations cannot capture the re-localization of the one-particle density matrix or the opposite-spin correlation. The RHF density matrix exhibits spurious delocalization in the stretched bond, and the RHF energy is well above the exact dissociation limit. Breaking symmetry in unrestricted Hartree-Fock (UHF) theory localizes the ↑- and ↓-spin electrons to different atoms, giving accurate dissociation energies at the expense of qualitatively incorrect spin densities. (Section SI-X of the supplementary material17 derives these effects for minimal basis H2 .) While H2 can be treated essentially exactly, strong correlation in larger systems remains a central challenge in electronic structure theory.68 Figure 8 shows that the EDR distinguishes the RHF, UHF, and FCI wavefunctions of stretched H2 . The figure presents contour plots of EDR(x; u) for the ↑-spin electron in stretched H2 . Here x is the position along the H–H bond axis, x = 0 is the bond midpoint, and u is the delocalization range. The RHF EDR exhibits unphysical delocalization around the bond midpoint, giving a peak near the bond midpoint that is shifted significantly upward to long delocalization length scales u. This is essentially identical to the EDR of stretched H+ 2 , which is shown in Figure SI-23 of the supplementary material.17 The UHF EDR shows the ↑-spin electron localized to the left H atom. This EDR is essentially identical to that of isolated H atom, shown in Figure SI-23 of the supplementary material.17 Figure 8(c) shows that the accurate FCI wavefunction puts 1/2 of the ↑-spin electron density on each atom, and gives a small EDR in the bonding region. This illustrates how the accurate wavefunction avoids both the spurious delocalization of RHF theory and the symmetry breaking of UHF theory. The small maximum at x = 0, u ∼ 5 bohrs arises in part from the ρ −1/2 (r ) prefactor in Eq. (3). The nonbonding character of the FCI wavefunction makes the FCI ρ(r ) much smaller than the corresponding RHF value for points r near the bond midpoint (see Figure SI-24 of the supplementary material17 ). For large u, the integrand of Eq. (3) remains modestly large at those points, leading to a peak from the ρ −1/2 (r ) contribution. Figure SI-24 of the


144104-8


J. Chem. Phys. 141, 144104 (2014)

FIG. 7. EDR(x, u) ((a), (c), (e)) and ELF ((b), (d), (f)) along the axes of the cyclohexatriene C=C bond ((a) and (b)), benzene C–C bond ((c) and (d)) and cyclohexatriene C–C bond ((e) and (f)). Gridlines are included to improve readability. Wavefunctions are from PBE0/6-31G(d,p) DFT, other details in text.

supplementary material17 also illustrates the EDR at two representative u, showing that the RHF EDR remains significantly larger than the FCI value. Plots of EDR(FCI, RHF; u) and EDR(FCI, UHF; u) (see Figure SI-25 of the supplementary material17 ) confirm that the UHF wavefunction is excessively localized, while the RHF wavefunction is excessively delocalized.

One notable result in Figure 8 is that the FCI EDR on each atom of stretched H2 is only ∼70% of the corresponding value on isolated H atom, such that the two peaks in Figure 8(c) are darker than the peak in Figure 8(b). This shows that the EDR detects the “fractional spin” states central to modern models of strong correlation.36, 72 Implications for using the EDR in such models are discussed below. Overall, we


144104-9


J. Chem. Phys. 141, 144104 (2014)

FIG. 8. EDR(x, u) of ↑-spin electrons in stretched H2 , evaluated from RHF (a), UHF (b), and nearly exact FCI (c) wavefunctions in the aug-cc-pVQZ basis.

suggest that this exploration of strong correlation in a model system illustrates the conceptual utility of the EDR. F. Excited states

Delocalization is central to understanding the properties of electronic excited states. In many-electron systems, excitations typically involve relatively small changes in the total wavefunction, such that analyses often consider differences between ground- and excited-state density matrices.73 We thus begin with a simple system in which excitation yields a large change in the density matrix. Figure 9 compares EDR(x; u) for the ground and first excited state of H2 at its equilibrium bond length. The figure follows the convention of Figure 8. As expected, the bonding→antibonding excitation reduces the amount of localized electrons at the center of the molecule, and increases the amount of delocalized electrons on the periphery of the molecule. Figure SI-26 of the supplementary material17 shows that EDR(ground, excited; u) is positive at all length scales u, indicating a small overlap between the excited-state density matrix and the test function in Eq. (2). We also consider the push-pull chromophore paranitroaniline (pNA). Its lowest intramolecular charge-transfer (CT) excitation is a π − π * transition involving charge transfer from NH2 to NO2 .74 Figure SI-27 of the supplementary material17 shows the calculated geometry. Figure 10 shows

the ground-state ELF and EDR, as well as the difference between ground and first CT excited states. We evaluate the excited-state ELF directly from the excited-state one-particle density matrix. As shown in Figure SI-19 of the supplementary material,17 these values are plotted in the π system, taken to lie in a plane 0.5 bohr above the nuclei.15 Figure 10(b) shows that the excited state ELF has a decrease in localized electrons above the NH2 nitrogen, and a substantial increase in localized electrons above the NO2 atoms. This is consistent with transfer of π electron density from NH2 to NO2 . Figure 10 shows that the EDR provides additional information about the pNA excitation. In the ground state, EDR(r ; u = 0.5 bohr) is significant in the NO2 group (Figure 10(c)), indicating that the NO2 π electrons are relatively localized. The difference between the ground- and excitedstate EDR(r ; u = 0.5 bohr) shows bright spots over the NO2 oxygens at u = 0.5 bohr (Figure 10(d)), indicating that the transferred π electrons are relatively localized. The difference between the ground- and excited-state EDR(r ; u = 0.5 bohr) shows dark spots over the NH2 nitrogen, and that of the EDR(r ; u = 1.0 bohr) shows dark spots over the NH2 hydrogens. This suggests that the electron density transferred from the NH2 is localized to ∼1 bohr. The difference between the ground- and excited-state ELF and EDR also show a dark spot over the carbon atom bound to NO2 . This is consistent with an orbital picture of the excitation: Figure 2 of Ref. 74

FIG. 9. EDR(x, u) in the ground state (a) and first excited state (b) of H2 . Wavefunction from RHF/6-311++G(2d,2p) calculations.


144104-10


J. Chem. Phys. 141, 144104 (2014)

FIG. 10. ELF ((a) and (b)), EDR(r ; u = 0.5 bohr) ((c) and (d)), and EDR(r ; u = 1.0 bohr)u = 1.0 bohr ((e) and (f)) for points r in the para-nitroaniline π system. Ground state GS ((a), (c), (e)) and difference between first CT excited and ground states ES-GS ((b), (d), (f)). Nitro and amino groups are at the left and right of the plot. Wavefunction from ωB97X/6-31+G(d,p) DFT calculations.

shows that this carbon atom has significant π -electron density in the HOMO, but has a node in the LUMO. The EDR difference plots in Figure 10 are still rather complicated, and more work is needed to refine their interpretation. However, they illustrate the EDR’s potential for interpreting the chemistry of excited states.

G. Delocalization and chemical reactivity

We close by showing how the EDR can build insight into chemical reactivity. Figure 11 illustrates the ELF and EDR in

the reactant, transition state, and product of the SN 2 reaction F − CH3 + Cl− → [F · · · CH3 · · · Cl]− → F− + H3 C − Cl. (5) The figure follows the convention of Figure 8, plotting EDR(x; u) where x is the position along the F–C–Cl bond axis and u is the delocalization length scale. The F–C–Cl bond is linear in all three structures. Figure SI-26 of the supplementary material17 shows the calculated geometries. The EDR plots in Figure 11 show a wealth of chemical information. The F and C core electrons, and the Cl n = 2


144104-11


J. Chem. Phys. 141, 144104 (2014)

FIG. 11. EDR(x, u) ((a)–(c)) and ELF ((d)–(f)) in the reactant, transition state, and product of Eq. (5). Results are plotted at positions x along the F–C–Cl bond axis. The fluorine nucleus is placed at the origin, the chlorine nucleus is near the right of the plot, and the carbon nucleus is marked with a vertical line in (a)–(c). The ordinate of (d)–(f) is the electron delocalization length scale u. Wavefunctions are from PBE0/6-311G(2d,2p) DFT, other details in text.

shell electrons, are visible at small u. The weak C–Cl− interaction on the right-hand side of the reactant (Figure 11(a)) is shifted downward to a more localized C–Cl bond in the late transition state (Figure 11(b)) and the product (Figure 11(c)). At the same time, the F–C bond on the left-hand side of the reactant (Figure 11(a)) bifurcates into a localized fluorine atom lone pair and delocalized electron density in the product H3 C region. Comparing the F–C bond at the left of Figure 11(a) to the C–Cl bond at the right of Figure 11(c) shows that the latter is shifted upward to longer delocalization lengths. This is consistent with increased delocalization of the weaker, longer bond to a second-row halogen. It is also consistent with the computed reaction energetics. Replacing the localized F–C bond with a more delocalized C– Cl bond is a thermodynamically unfavorable reaction with a positive reaction energy (26.1 kcal mol−1 , Ref. 75) and a late transition state. The ELF shows that the reactant’s F–C bond (Figure 11(d)) bifurcates into two peaks in the transition state (Figure 11(e)), while the C–Cl interaction coalesces into a single peak (Figure 11(f)). However, the ELF does not directly quantify the relative delocalization of the forming and breaking bonds. This treatment of a simple model system motivates application of the EDR to more complicated and less wellcharacterized reactions. V. DISCUSSION

The electron delocalization range complements existing electronic structure descriptors1, 2 by visualizing and quantifying the delocalization of electrons in different atomic shells (Figure 2(a)), valence electrons of different atoms (Figure 3),

solvated and weakly bound electrons (Figure 4), and electrons in metals and insulators (Figure 6(c)). The EDR highlights the different degrees of delocalization present in solvent cavities (Figure 5), in stretched or compressed bonds (Figure 7), in excited states (Figures 9 and 10), and in transition states (Figure 11). It shows potential for visualizing the interplay of delocalization and strong correlation (Figure 8). The EDR does all of this at modest computational cost. The EDR is built from inexpensive one-electron operators. It requires only overlap integrals between atomic orbitals and gridpointcentered auxiliary basis functions, integrals which can be evaluated analytically using standard algorithms.16, 20 This work lays a foundation for applying the EDR across the chemical enterprise. We envision using the EDR to visualize delocalization and strong correlation in electrides,56–58 metal-insulator phase transitions,76 and the ground and excited states of defects in solids and on surfaces. To illustrate, Figure 12 shows how the EDR distinguishes a chlorine atom vacancy in solid NaCl.77 (Figure SI-27 of the supplementary material17 shows the calculated unit cell.) We will use the EDR to quantify the extent of electron localization and the role of d electrons in transition metal catalysis, as well as the delocalization of excited states.73 Evaluating the gradient field ∇EDR(r ; u) at different length scales u will build on the enormously successful topological analyses of the electron density and ELF.9, 11 Replacing γ (r , r ) in Eq. (2) with its particle number derivative (the Fukui matrix78 ) offers potential new insights into the delocalization of molecules’ most reactive electrons. The EDR’s detection of fractional occupancy illustrated in Figure 8 also motivates applications to real-space models of strong correlation.36–38


144104-12


J. Chem. Phys. 141, 144104 (2014)

FIG. 12. ELF (a), EDR(r ; u = 1 bohr) (b), and EDR (r ; u = 4 bohrs) (c) in bulk NaCl. Results are plotted in a plane containing a Cl atom vacancy: Cl atom is at lower left, Na atoms are at upper left and lower right, vacancy is at upper right. Wavefunction is from LSDA/STO-3G DFT. The high defect concentration (1 out of every 4 Cl atoms removed) suffices for illustration.

We also envision technical improvements. Replacing gu (r , r ) in Eq. (2) with test functions such as exp(−(|r − r |−u)2 ) should provide a sharper probe of specific length scales u. This may ameliorate the “smearing” of the different shells visible in Figure 2, and the loss of structure in the valence region visible in Figure 4. Non-spherically-symmetric test functions characterizing electron delocalization along different directions could provide insights into anisotropic systems such as graphene nanoribbons or conjugated polymers. All of these test functions may be expanded in auxiliary basis sets as described in Ref. 16. In conclusion, the EDR’s answer to the question “How far do electrons delocalize?” shows promise for quantifying, visualizing, and building insight into the fundamentally nonclassical aspects of chemical bonding and electronic structure.

ACKNOWLEDGMENTS

B.G.J. acknowledges support by the Qatar National Research Foundation through the National Priorities Research Program (NPRP Grant No. 09-143-1-022), and by the Department of Chemistry at Texas Christian University. The authors thank Marcel Nooijen for comments on strongly correlated systems. 1 J.

Poater, M. Duran, M. Solá, and B. Silvi, Chem. Rev. 105, 3911 (2005). Matito and M. Solá, Coord. Chem. Rev. 253, 647 (2009). 3 G. N. Lewis, Valence and the Structure of Atoms and Molecules (Dover, New York, 1966). 4 A. E. Reed, R. B. Weinstock, and F. Weinhold, J. Chem. Phys. 83, 735 (1985). 5 W. Yang and R. G. Parr, Proc. Natl. Acad. Sci. U.S.A. 82, 6723 (1985). 6 P. Pérez, E. Chamorro, and P. W. Ayers, J. Chem. Phys. 128, 204108 (2008). 7 A. D. Becke and K. E. Edgecombe, J. Chem. Phys. 92, 5397 (1990). 8 A. Savin, O. Jepsen, J. Flad, O. K. Andersen, H. Preuss, and H. G. von Schnering, Angew. Chem. Int. Ed. Engl. 31, 187 (1992). 9 B. Silvi and A. Savin, Nature (London) 371, 683 (1994). 10 A. Savin, J. Mol. Struct. 727, 127 (2005). 11 R. F. W. Bader, Atoms in Molecules: A Quantum Theory (Oxford University Press, Oxford, 1990). 12 Spin and time dependence are generally suppressed for conciseness. All wavefunctions, densities, density matrices, and so on are for majority-spin electrons in stationary states of Schrödinger’s equation, or in the linear response regime of an infinitesimal time-dependent perturbation. 2 E.

13 R.

S. Mulliken, J. Chem. Phys. 23, 1833 (1955). Levy, Phys. Rev. A 43, 4637 (1991). 15 R. F. W. Bader, A. Streitwieser, A. Neuhaus, K. E. Laidig, and P. Speers, J. Am. Chem. Soc. 118, 4959 (1996). 16 B. G. Janesko, G. Scalmani, and M. J. Frisch, J. Chem. Phys. 141, 034103 (2014). 17 See supplementary material at http://dx.doi.org/10.1063/1.4897264 for implementation details, derivations for stretched H2 , and Figures SI-1– SI-27. 18 R. Krishnan, J. Binkley, R. Seeger, and J. Pople, J. Chem. Phys. 72, 650 (1980). 19 T. H. Dunning, Jr., J. Chem. Phys. 90, 1007 (1989). 20 P. M. W. Gill, in Advances in Quantum Chemistry (Academic Press, San Diego, CA, 1994), Vol. 25, pp. 141–205. 21 Distinctions among, e.g., “core” and “valence” electrons are a colloquialism used here to aid understanding. All calculations treat electrons as indistinguishable Fermions. 22 E. R. Johnson, S. Keinan, P. Mori-Sánchez, P. Contreras-García, A. J. Cohen, and W. Yang, J. Am. Chem. Soc. 132, 6498 (2010). 23 K. Wagner and M. Kohout, Theor. Chem. Acc. 128, 39 (2011). 24 H. L. Schmider and A. D. Becke, J. Mol. Struct.: THEOCHEM 527, 51 (2000). 25 M. Kohout, Int. J. Quantum Chem. 97, 651 (2004). 26 R. Ponec and J. Roithová, Theor. Chem. Acc. 105, 383 (2001). 27 R. Ponec and D. L. Cooper, J. Phys. Chem. A 111, 11294 (2007). 28 R. Ponec, D. L. Cooper, and A. Savin, Chem. Eur. J. 14, 3338 (2008). 29 R. F. W. Bader and M. E. Stephens, J. Am. Chem. Soc. 97, 7391 (1975). 30 J. Cioslowski and S. T. Mixon, J. Am. Chem. Soc. 113, 4142 (1991). 31 X. Fradera, M. A. Austen, and R. W. F. Bader, J. Phys. Chem. A 103, 304 (1999). 32 B. Silvi, Phys. Chem. Chem. Phys. 6, 256 (2004). 33 H. L. Schmider, J. Chem. Phys. 105, 11134 (1996). 34 H. L. Schmider and A. D. Becke, J. Chem. Phys. 116, 3184 (2002). 35 E. Proynov, F. Liu, and J. Kong, Phys. Rev. A 88, 032510 (2013). 36 A. D. Becke, J. Chem. Phys. 119, 2972 (2003). 37 A. D. Becke, J. Chem. Phys. 122, 064101 (2005). 38 A. D. Becke, J. Chem. Phys. 138, 074109 (2013). 39 A. D. Becke and M. R. Roussel, Phys. Rev. A 39, 3761 (1989). 40 B. G. Janesko and G. E. Scuseria, J. Chem. Phys. 127, 164117 (2007). 41 B. G. Janesko and G. E. Scuseria, J. Chem. Phys. 128, 084111 (2008). 42 M. J. Frisch, G. W. Trucks, H. B. Schlegel et al., Gaussian Development Version, Revision H.35, Gaussian, Inc., Wallingford, CT, 2010. 43 F. E. Jorge, E. V. R. de Castro, and A. B. F. da Silva, Chem. Phys. 216, 317 (1997). 44 R. G. Parr and W. Yang, Density-Functional Theory of Atoms and Molecules (Oxford University Press, New York, 1989). 45 W. Kohn and L. Sham, Phys. Rev. 140, A1133 (1965). 46 A. Seidl, A. Görling, P. Vogl, J. A. Majewski, and M. Levy, Phys. Rev. B 53, 3764 (1996). 47 S. H. Vosko, L. Wilk, and M. Nusair, Can. J. Phys. 58, 1200 (1980). 48 J. P. Perdew, K. Burke, and M. Ernzerhof, Phys. Rev. Lett. 77, 3865 (1996). 49 C. Adamo and V. Barone, J. Chem. Phys. 110, 6158 (1999). 14 M.


144104-13 50 M.


Ernzerhof and G. E. Scuseria, J. Chem. Phys. 110, 5029 (1999). Chai and M. Head-Gordon, J. Chem. Phys. 128, 084106 (2008). 52 Specialists will recall that generalized Kohn-Sham DFT actually provides the wavefunction of a reference system of N noninteracting Fermions, moving in a nonlocal one-electron potential, with the same density as the real system.45, 46 53 N. C. Handy and H. F. Schaefer III, J. Chem. Phys. 81, 5031 (1984). 54 V. I. Lebedev, Zh. Vychisl. Mat. Mat. Fiz. 16, 293 (1976). 55 L. Turi and P. Rossky, Chem. Rev. 112, 5641 (2012). 56 J. M. Dye, Acc. Chem. Res. 42, 1564 (2009). 57 E. R. Johnson, A. Otero-de-la Roza, and S. G. Dale, J. Chem. Phys. 139, 184116 (2013). 58 S. G. Dale, A. Otero-de-la Roza, and E. R. Johnson, Phys. Chem. Chem. Phys. 16, 14584 (2014). 59 L. Kevan, Acc. Chem. Res. 14, 138 (1981). 60 O. A. Vydrov and G. E. Scuseria, J. Chem. Phys. 125, 234109 (2006). 61 O. A. Vydrov, J. Heyd, A. V. Krukau, and G. E. Scuseria, J. Chem. Phys. 125, 074106 (2006). 62 At this level of theory, graphene and boron nitride are respectively predicted to have indirect bandgaps 0.28 and 4.51 eV. This imperfect treatment of the graphene band structure suffices for illustrative purposes. 63 S. N. Steinmann, D. F. Jana, J. I.-C. Wu, P. von R. Schleyer, Y. Mo, and C. Corminboeuf, Angew. Chem., Int. Ed. 48, 9828 (2009). 51 J.-D.

J. Chem. Phys. 141, 144104 (2014) 64 J. Zhu, C. Dahlstrand, J. R. Smith, S. Villaume, and H. Ottosson, Symmetry

2, 1653 (2010). Kalescky, E. Kraka, and D. Cremer, J. Phys. Chem. A 118, 223 (2014). 66 K. B. Wiberg, Acc. Chem. Res. 29, 229 (1996). 67 J. P. Perdew, A. Ruzsinszky, L. A. Constantin, J. Sun, and G. I. Csonka, J. Chem. Theory Comput. 5, 902 (2009). 68 K. Burke, J. Chem. Phys. 136, 150901 (2012). 69 A. D. Becke, J. Chem. Phys. 140, 18A301 (2014). 70 O. Gunnarsson and B. I. Lundqvist, Phys. Rev. B 13, 4274 (1976). 71 D. Cremer, Mol. Phys. 99, 1899 (2001). 72 A. J. Cohen, P. Mori-Sánchez, and W. Yang, J. Chem. Phys. 129, 121104 (2008). 73 S. Tretiak and S. Mukamel, Chem. Rev. 102, 3171 (2002). 74 J. J. Eriksen, S. P. A. Sauer, K. V. Mikkelsen, O. Christiansen, H. J. A. Jensen, and J. Kongsted, Mol. Phys. 111, 1235 (2013). 75 Y. Zhao, B. J. Lynch, and D. G. Truhlar, J. Phys. Chem. A 108, 2715 (2004). 76 S. Sharma, J. K. Dewhurst, S. Shallcross, and E. K. U. Gross, Phys. Rev. Lett. 110, 116403 (2013). 77 G. Mallia, R. Orlando, C. Roetti, P. Ugliengo, and R. Dovesi, Phys. Rev. B 63, 235102 (2001). 78 P. Bultinck, D. Clarisse, P. W. Ayers, and R. Carbo-Dorca, Phys. Chem. Chem. Phys. 13, 6110 (2011). 65 R.


Twenty years of Milan criteria: how far do we go.

Targeted therapies for pancreatic adenocarcinoma: Where do we stand, how far can we go?

How far can mice carry autism research?

Ebola Vaccine: How Far are we?

2014-2016: How far has LUGPA come?

How far to lower blood pressure.

Thyroid cytopathology: how far can we go?

How far should blood pressure be lowered?

Editorial: How do you do, Mrs. Jones?

How I Do It.

How do oocytes disappear?

How do I look?

How do phagocytes eat?

How do you feel?

Prevention of diabetes: How far have we gone?

Realistically, how far are we from a universal malaria drug?

How far from the SNP may the causative genes be?

Active Surveillance in Prostate Cancer: How Far Should We Go?

Simplifying atrial fibrillation ablation: how far can we go?

Pharmacogenetics of smoking: how far to the clinic?

Retrospective studies and pancreatic adenocarcinoma: how far can we backdate?

Deaths from falls in children: how far is fatal?

[Treatment of early onset scoliosis : How far can we go?].

Toward Precision Medicine: How Far Is the Goal?