Publication Cover
Molecular Physics
An International Journal at the Interface Between Chemistry and Physics
Volume 114, 2016 - Issue 7-8: Special Issue in honour of Andreas Savin
425
Views
5
CrossRef citations to date
0
Altmetric
Development and Application of Quantum-Chemistry Interpretative Methods

Atomic shell structure determined by the curvature of the electron position uncertainty

Pages 1297-1303 | Received 30 Sep 2015, Accepted 25 Oct 2015, Published online: 24 Nov 2015

abstract

The electron density logarithm is closely connected to the Shannon information entropy characterising the spread of the electron density. It can be seen as the uncertainty to predict the position of an electron. The position uncertainty curvature, defined as the negative of the Laplacian of the electron density logarithm, describes regions of space where the uncertainty of the electron position is higher, respectively lower, than the average in the surrounding. The atomic shells can be located in regions of high-position uncertainty curvature, thus resolving the atomic shell structure of the atoms Li to Xe. The shell boundaries are given as minima of the uncertainty curvature. This indicator, based alone on the electron density, is suitable to describe the bonding situation in molecules and solids.

GRAPHICAL ABSTRACT

1. Introduction

The spherically averaged electron density ρ(r) of an atom displays a maximum at the atomic position followed by monotonically exponential decay [Citation1–5]. The idea of atomic shells successively filled up by electrons, as anticipated from the Periodic Table of elements, is not directly reflected by this simple distribution. To accomplish this, specific functions derived from the electron density or even functions based on the orbital representation of wave function are needed. The analysis of such functions revealed, at least partially, the atomic shell structure [Citation2,Citation4,Citation6–22]. However, only few functionals were able not only to resolve the atomic shell structure, but also to yield reasonable electron populations for the shell occupations, for instance the average local electrostatic potential [Citation18], the electron localisation function (ELF) [Citation14,Citation17], as well as ELF modified by the incorporation of the electrostatic potential [Citation19], respectively, the electron localisability indicator (ELI) [Citation23–25] and the one-electron potential (OEP) [Citation20]. The electron localisability indicator variant ELI-D [Citation26] yields at the Hartree–Fock (HF) level in principle the same topology as ELF, but in contrast to ELF is defined also at correlated level of theory [Citation26,Citation27].

The piecewise exponential decay of the electron density ρ(r) of an atom was investigated by Wang and Parr [Citation4] who showed that for the ground state of atoms, the function logρ(r) reflects shell structure by a significant change of the logρ(r) slope. This circumstance was used later for the quantity -|ρ(r)|/ρ(r) displaying separate plateaus for each atomic shell [Citation16]. For the quantitative evaluation of the atomic shells, suitable separators were needed to determine the shell boundaries of -|ρ(r)|/ρ(r). Bohórquez and Boyd, after connecting -ρ(r)/ρ(r) with the quantum spread, used the radial inflection points as the shell delimiters [Citation21]. Of course, such definition based on spherical symmetry is in principle applicable only to shells of free atoms. In the present study, a function derived from the electron density logarithm was used to visualise and analyse the shell structure of atoms. The function can easily be applied also to the analysis of bonding situation in molecules and solids.

2. Theory

For the quantity X of n discrete values xi measured with the respective probabilities pi the expression: (1) Sp=-inpilogpi(1) is the Boltzmann–Shannon information entropy of the measurement [Citation28]. The information entropy describes the uncertainty to obtain one (further not specified) discrete value of the quantity X. If specific event has the probability pk = 1 and the remaining ones equal zero, then the information entropy Sp = 0. Of course, in this case, the result of the measurement will be certainly xk. On the other hand, if all n values occur with the same probability pi = 1/n the uncertainty about the outcome of the measurement reaches its maximum and the information entropy Sp = log n. The information entropy as given by Equation (Equation1) can also be seen as the average value of the quantity Y of n discrete values yi = −log pi occurring with the probabilities pi. Each single value −log pi can be regarded as the uncertainty to obtain the value xi for the quantity X.

The electron density ρ(r) is connected with the probability (multiplied by the number of electrons in the system) of finding the coordinates of an electron between r and dr, with the remaining electrons distributed anywhere in the space. Let us take for the measured quantity the determination of the position r of an electron. The information entropy Sρ for this measurements, describing the spread of the electron density, is given by: (2) Sρ=-ρ(r)logρ(r)dr.(2) In accordance with the discrete case, the function -logρ(r) can be interpreted as the uncertainty in the prediction of the electron position and the information entropy Sρ as the mean value of this uncertainty [Citation29]. Often the terms local Shannon information for s(r)=-ρ(r)logρ(r) as well as Shannon information density per particle for the function -logρ(r) are used [Citation30]. Shannon entropy was also connected with the correlation strength [Citation31] as well as utilised in the orbital-free density-functional theory (DFT) [Citation32]. The combined Shannon entropies in position and momentum space are involved in the derivation of the Heisenberg uncertainty relation [Citation33].

The investigation of Wang and Parr [Citation4] and Kohout et al. [Citation16] corroborated the idea that the electron position uncertainty -logρ(r) as well as its gradient -ρ(r)/ρ(r), also termed the local wave-vector, is able to detect the shell structure of atoms. However, as Wang and Parr stated ‘... the transitions from one exponential to the next occur over certain intervals and cannot be associated with single points.’ Some years later, Bohórquez and Boyd connected -ρ(r)/ρ(r) with the so-called quantum spread P˜, which is a part of the local representation of the momentum operator [Citation21]: (3) P˜(r)=-2ρ(r)ρ(r).(3) Analysing the representation of the atomic shell structure with the magnitude of the quantum spread P˜, the above authors defined the radial inflection points 2P˜/r2=0 located after each local maximum as the shell boundaries. The application of this procedure for some closed shell atoms yielded shell populations close to the one of ELF, with somewhat lower populations of the first shell when compared to the ELF results (which overestimate the first shell population) [Citation21]. This approach could easily be applied also to atoms with spherically averaged electron density. However, it would be more complex to employ such definition of the shell boundaries for molecules due to the los of the sphericity.

This obstacle can be simply circumvented by the utilisation of the Laplacian of the electron position uncertainty, instead of the radial curvature. Let us define the position uncertainty curvature (PUC): (4) PUC =-2logρ(r)=ρ(r)ρ(r)2-2ρ(r)ρ(r).(4) The PUC is given by two term – the scaled density gradient and the scaled density Laplacian, respectively – in ratio 1:1. Interestingly, there is another density function also composed of these two terms, namely the OEP defined by 1/2(2ρ)/ρ [Citation12,Citation34]. For OEP, the corresponding ratio of the two terms is 1:2. The OEP was used to visualise the chemical bonding [Citation12] and it nicely resolves the atomic shell structure [Citation20].

In case of spherical symmetry, the ∇2 operator is replaced by ∂2/∂r2 + (2/r)∂/∂r, i.e. the second right-side term of Equation (Equation4) will be split into two contributions. For large distances from the nucleus ρ(r) ≈ krβe−αr (with k, α > 0 and β ≤ 0) [Citation35]. Applying this to Equation (Equation4) yields 2α/r − β/r2, i.e. the PUC asymptotically approaches zero from above for infinite distance r. At the origin ,∇2ρ(r) is negative infinite, whereas the first term on the right hand side of Equation (Equation4) approaches a constant (see for instance Ref. [Citation10]). Hence, a positive infinite PUC is found at the atomic origin.

A positive curvature is found at such points, where the position uncertainty is lower than the average value in the closest neighbourhood. The uncertainty in the localisation of an electron (compared to the average value in the closest neighbourhood) decreases, i.e. the electron density is less spread, with increasing value of PUC. Thus, it is an appealing idea to assign the (positive) peaks in a -2logρ(r) vs. radius diagram to the atomic shells and the positions of the minima between the peaks to the shell radii. This definition of shell radii based on PUC is reasonable also from another point of view. Namely, because for the second row atoms, two shells are expected (and two PUC maxima are indeed found, cf. (a) for the Be atom). Whereas the first shell has a finite radius (correspondingly, there is a negative PUC minimum), the second one extends to infinity. In contrast, choosing the zeros of PUC as the shell delimiters would result in two finite shell radii (otherwise one of the odd or even zeros had to be chosen), cf. (a). Because PUC is defined independently of the symmetry of the systems, it can be easily applied to molecules and crystals.

Figure 1. Atomic shells represented by PUC using the wave functions of Clementi and Roetti. (a) Be atom. (b) Ag atom. The spurious minimum at 4.43 bohr within the 5th shell is marked by the arrow.

Figure 1. Atomic shells represented by PUC using the wave functions of Clementi and Roetti. (a) Be atom. (b) Ag atom. The spurious minimum at 4.43 bohr within the 5th shell is marked by the arrow.

3. Results and discussion

Using the wave functions of Clementi and Roetti [Citation36] for the atoms Li to Xe, the spherically averaged electron densities ρ(r) as well as the PUC = −∇2log ρ were calculated with the DGrid program [Citation37]. The minima of PUC were detected and assigned to the atomic shell radii. The electron density within the radial shells were integrated yielding the corresponding shell populations. The results are compiled in .

Table 1. The shell radii and electron populations given by the minima of −∇2log (ρ).

The shell radii, and thus also the shell populations, are very close to the data obtained for ELF using the so-called ‘closed-shell’ formula [Citation17] (at the single-determinantal level this applies also to ELI-D). The shell populations of the first three shells as given by PUC differ at most by 0.1 electrons from the data determined by the ELF. Correspondingly, this means that also for PUC the first shell is slightly over-populated by 2.3 electrons, starting already at the V atom. The same is true for the second shell occupied up to 8.8 electrons in case of Cd to Xe. In contrast, up to 1 electron is missing in the third shell (in case of full occupation according to the Periodic Table). In case of PUC, the shell occupations of the fourth shell are closer to 18 electrons (for full shell occupation) than are the corresponding shell populations determined for the ELF.

There are four exceptions where the valence shells are not properly resolved. In case of the Pd atom, the fifth shell is not occupied according to the Aufbau principle. However, PUC shows an unexpected minimum for the shell boundary with 0.7 electrons present in this spurious ‘valence’ shell. Indeed, this behaviour is similar to the one for atomic shells determined from the OEP, in which case even 2.3 electrons are found in such fifth shell for the Pd atom [Citation20]. The another exception is the Ag atom for which the ‘ideal’ shell occupation [Citation38] assumes 1 electron in the valence shell. Similarly, ELF delivers 1.1 electron for the fifth shell, whereas the OEP valence shell is over-populated by 2.8 electrons [Citation20]. Although there are 0.9 electrons outside the fourth PUC shell for the Ag atom, cf. , the valence shell is split into two parts. As can be seen from (b), there is a shallow PUC minimum at 4.43 bohr which formally splits the valence shell of Ag into a part populated by 0.6 electrons and an outer part populated by 0.3 electrons.

The two remaining exceptions concern a missing valence shell. For the atoms Ru and Rh, the Aufbau principle predicts valence shells occupied by 1 electron. The ELF resolves the valence shells for both atoms, with reasonable occupations of 1.3 and 1.2 electrons for Ru and Rh, respectively. Even the OEP shows separate valence shell for Ru and Rh, although heavily over-populated by almost 3 electrons. In contrast, there are no populations given in for the fifth shell of Ru and Rh. The inspection of helps to clarify this issue. Whereas all four inner shells of Ru and Rh are well represented by PUC maxima, there is in both cases only a shoulder indirectly suggesting the presence of the valence shell. Not even the utilisation of a wave function from, for instance, a relativistic density functional calculation will amend this feature.

Figure 2. Atomic shells represented by PUC using the wave functions of Clementi and Roetti. The arrows mark the shoulder present instead of the 5th shell maximum. (a) Ru atom. (b) Rh atom.

Figure 2. Atomic shells represented by PUC using the wave functions of Clementi and Roetti. The arrows mark the shoulder present instead of the 5th shell maximum. (a) Ru atom. (b) Rh atom.

Proper representation of the atomic shell structure is an essential prerequisite for plausible analysis of the bonding situation in molecules and solids. Then, the core regions which do not strongly participate on the bonding can be reasonably separated in the position space from the valence region that can be subsequently analysed in detail. As an example, the bonding of the cyclopropane C3H6 was analysed with PUC. The C3H6 was computed with the program GAMESS [Citation39] at the HF level using the correlation consistent triplet zeta basis set CCT. The optimised HF distances amount to d(CC) = 2.828 bohr and d(CH) = 2.028 bohr, respectively. With the program DGrid, the output of GAMESS was converted into a wave function data file and the PUC distribution computed on an equidistant grid of 0.025 bohr mesh size.

In (a), the spherical isosurfaces of PUC = 4.25 marking the core regions of the carbon atoms can be observed. The corresponding core PUC basins enclose 2.09 electrons. This is a situation closely resembling the ELI-D results for the HF calculation of C3H6. The 4.25-localisation domains include the almost spherical isosurfaces around the hydrogen positions. Additionally, there are also isosurfaces located between the C and H isosurfaces along the C–H line. Each hydrogen PUC basin encloses 1.32 electrons, whereas the C–H basin is populated by 0.75 electrons. This is a situation very different from the one for ELI-D (respectively, ELF), which typically shows only a single hydrogen localisation domain with an ELI-D basin populated (in case of C3H6) by 2.03 electrons. A separate ELI-D domain between the C and H atoms is not present. Similarly, there are two localisation domains, displayed in (a) as red coloured isosurfaces of PUC = 3.0, between the carbon atoms. In case of ELI-D, there is only single localisation domain (single attractor) representing the C–C bond. This behaviour can be easily explained by the inspection of the PUC and ELI-D topology for the free atoms. Whereas the positions of the minima (the atomic shell boundaries) are almost identical, the positions of the shell maxima for PUC (1.028 bohr) and ELI-D (1.488 bohr) are very different for the C atom. Taking as a first approximation for the C–C bond a simple overlap of the atomic PUC and ELI-D distributions shows that the PUC maxima are much closer to the nuclei than the C–C midpoint (1.414 bohr), whereas the ELI-D maxima are located slightly behind the C–C midpoint. This is a typical topological feature of PUC showing separate shell maxima and corresponding saddle point as the signature of bond. Interestingly, this allows in certain sense to define kind of bond polarity, because the bond is represented by two PUC basins related to separate participating atom. Of course, in case of the C–C bond, each of the C atoms yields the same contribution (0.88e) to the population of the bond (super)basin. Similarly, considering the combined hydrogen and CH PUC-basin (together populated by 2.07e) as a ‘split’ ELI-D hydrogen basin would yield 64% of the combined population to be located at the H atom and 36% attributed to the C atom. This seems to be far too polar for the C–H bond in the cyclopropane.

Figure 3. PUC localisation domains for the C3H6 and N2 molecules. (a) HF calculation of C3H6 using the CCT basis. Red coloured: isosurfaces of PUC = 3.0. Yellow coloured: isosurfaces of PUC = 4.25. (b) HF calculation of N2 using the CCT basis. Blue coloured: 1.77-localisation domain. Yellow coloured: isosurfaces of PUC=4.4. (c) CISD calculation of N2 using the CCT basis. Yellow coloured: 4.4-localisation domains. (d) DFT calculation of N2 using the TZP basis. Yellow coloured: 4.2-localisation domains.

Figure 3. PUC localisation domains for the C3H6 and N2 molecules. (a) HF calculation of C3H6 using the CCT basis. Red coloured: isosurfaces of PUC = 3.0. Yellow coloured: isosurfaces of PUC = 4.25. (b) HF calculation of N2 using the CCT basis. Blue coloured: 1.77-localisation domain. Yellow coloured: isosurfaces of PUC=4.4. (c) CISD calculation of N2 using the CCT basis. Yellow coloured: 4.4-localisation domains. (d) DFT calculation of N2 using the TZP basis. Yellow coloured: 4.2-localisation domains.

(b) shows the PUC for the N2 molecule computed with GAMESS at the HF level using the CCT basis set at the optimised distance of 2.017 bohr and evaluated on a grid of 0.025 bohr mesh size. Again, the core region can be nicely recognised by the two spherical 4.4-localisation domains. The correspond core basins are populated by 2.1 electrons. In the valence region, the two lone pair regions are visible by separate outwards pointing 4.4-localisation domains with the corresponding basins populated by 3.2 electrons. There is also a separate irreducible 4.4-localisation domain for the bond. Thus, for the HF wave function, there is only single PUC attractor at the bond line. The corresponding basin is populated by 2.9 electrons. Additionally, the blue coloured 1.77-localisation domain located around a ring attractor describes a second bond feature, with the basin populated by 0.4 electrons. Besides this ring attractor, the description of the nitrogen bond is similar to the one with ELI-D, also showing with 2.9e somewhat over-populated lone pair basins.

For the N2 molecule, a configuration-interaction calculation was performed with single and double excitations (CISD) of 10 electrons (i.e. with inactive core electrons) into 58 orbitals using the CCT basis set (the optimised bond distance of 2.059 bohr was used). The PUC distribution shown in (c) displays at first glance less features than the diagram for the HF calculation, because of the absence of the bonding ring localisation domain. Otherwise, the 4.4-localisation domains again mark the core regions and lone pair regions as well as the bonding region. The corresponding core, lone pair, and (combined) bond basins are populated by 2.1e, 3.8e, and 2.1e, respectively. This means that for the CISD wave function, the lone pair as represented by the PUC basin is even more over-populated than in case of the HF wave function. In contrast to the HF PUC topology, as already indicated above by the word ‘combined’, the bond basin for the CISD calculation is split into two basins, i.e. the 4.4-localisation domain in (c) is reducible due to a shallow PUC saddle point at the bond midpoint.

For comparison, the DFT calculation of N2 was performed with the ADF program [Citation40–42] using the TZP basis set and the BLYP functional at the bond distance of 2.060 bohr. The PUC was computed on a grid of 0.025 bohr mesh size. (d) with the 4.2-localisation domains shows the core and lone pair regions together with clearly separated localisation domains representing the bond. The populations of the corresponding basins is almost the same as the for the CISD calculation. The pronounced separation of the two bond PUC maxima results from changes in the electron density. The comparison between the electron density for the DFT and CISD calculation reveals slight decrease of the density in the bond and lone pair region accompanied by a density increase in the inter-shell regions. The density decrease around the bond midpoint changes the curvature of the electron density, thus, affecting the distribution of PUC.

4. Conclusion

The Shannon information entropy derived from the electron density describes the spread of the density. It can be seen as the mean value of the Shannon information density per particle -logρ(r), which is connected with the uncertainty in the prediction of the position of an electron. The position uncertainty can be inspected from the topological viewpoint by the analysis of its Laplacian. This defines the PUC = -2logρ(r), that indicates whether the uncertainty is higher or lower than the average in the surrounding. PUC is able to resolve the atomic shell structure for almost all atoms from Li to Xe not only qualitatively but also to yield reasonable shell populations, similar to the one for ELI-D/ELF. The indicator PUC is applicable to molecules and crystals. It was demonstrated that localisation domains of PUC can be used as bonding descriptors representing core region, as well as lone pairs and bonds. Because the PUC atomic shell maxima are located closer to the corresponding nuclei, the bonds are usually characterised by two irreducible localisation domains and a PUC saddle point. The PUC is derived from the electron density and thus can be evaluated also from the experimental densities.

Disclosure statement

No potential conflict of interest was reported by the author.

References