Theoretical study of atoms by the electronic kinetic energy density and stress tensor density Hiroo Nozaki, Kazuhide Ichikawa, and Akitomo Tachibana∗ Department of Micro Engineering, Kyoto University, Kyoto, 615-8540, Japan Abstract 6 1 0 We analyze the electronic structure of atoms in the first, second and third periods using the 2 n electronic kinetic energy density and stress tensor density, which are local quantities motivated a J by quantum field theoretic consideration, specifically the rigged quantum electrodynamics. We 2 2 compute the zero surfaces of the electronic kinetic energy density, which we call the electronic ] interfaces, of the atoms. We find that their sizes exhibit clear periodicity and are comparable h p to the conventional atomic and ionic radii. We also compute the electronic stress tensor density - m and its divergence, tension density, of the atoms, and discuss how their electronic structures are o t a characterized by them. . s c i s PACS numbers: 31.15.-p, 32.10.-f, 32.90.+a y h Keywords: atoms, kinetic energy density, stress tensor density, tension density, core-valence partitioning p [ 1 v 5 8 8 5 0 . 1 0 6 1 : v i X r a ∗Electronic address: [email protected] 1 I. INTRODUCTION In order to understand the quantum systems from the viewpoint of quantum field theory, in particular quantum electrodynamics (QED), one of the authors has developed the rigged QED (RQED) theory [1–5] and our group has applied it to various molecular and condensed matter systems [1–33]. Due to the field theoretic nature, the theory has local quantities defined at each point in space which are useful to describe quantum systems. Most im- portant quantities are the electronic kinetic energy density and stress tensor density, whose definitions and meanings are presented in the next section. In this section, we shall review our past studies using them. As for the electronic kinetic energy density, our definition is proposed in Ref. [1], which is motivated by the quantum field theoretic consideration, and one of its features is that it is not positive-definite. Then, it has been proposed that the outermost zero surface of the kinetic energy density gives the intrinsic shape of atoms and molecules. This zero surface is designated as “electronic interface” in Ref. [1], and a variety of chemical reactions has been studied using it [6–12, 14–16, 19, 20, 24, 28–32]. We briefly describe some of them below. In Ref. [15], the volume surrounded by the outermost electronic interfaces of a cluster model is defined as its “shape volume” and is used to calculate the electronic dielectric constant of the cluster models of silicate compounds by means of the Clausius-Mossotti equation. In Refs. [20] and [31], local reactivity of Pt clusters and Al clusters regarding hydrogen adsorption has been investigated by mapping the regional chemical potential [1, 17, 18, 34], which measures the local chemical reactivity, over the electronic interface of the clusters. The same technique is used to study adsorption of lithium atom on carbon nanotube in Ref. [30]. Meanwhile, using the electronic stress tensor density, together with its divergence, which iscalledtheelectronictensiondensity, newviewsonchemicalbondinghasbeenproposed[1– 5,11,13]. Wehaveshownthatthe“Lagrangepoint”, thepointbetweentwoatomswherethe tension density vanishes, can well characterize the chemical bond between them [17, 18]. We have also shown that the “Lagrange surface”, separatrices of the tension density vector field, may define boundary surfaces of atoms in a molecule [3, 4, 26, 27, 33]. It is pointed out that the eigenvalues and corresponding eigenvectors of the electronic stress tensor in the bonding region play important roles for characterizing types of chemical bonding such as covalency 2 [11, 13] and metallicity [4, 17, 32]. In addition, the bonds between some of semimetal atoms are found to have intermediate nature between covalent and metallic bonds in view of the stress tensor density [33]. Although we have applied the analyses based on the electronic kinetic energy density and stress tensor density to molecular and periodic systems as above, there have not been a systematic application to atoms. It goes without saying that atoms are building blocks for molecules and their electronic kinetic energy density and stress tensor density carry very important information for our analysis. Hence, in this paper, we compute the electronic kinetic energy density, stress tensor density, and tension density of the atoms in the first, second and third periods, and investigate how the atoms are characterized by them. This paper is organized as follows. In Sec. II, we show our definition of electronic kinetic energy density and stress tensor density. We also describe our computational setups. In Sec. IIIA, we investigate the shape of the electronic interfaces of atoms in the first, second and third periods. In Sec. IIIB, we study the size of the electronic interfaces and compare with the conventionally defined atomic and ionic radii. In Sec. IIIC, we study the integrated electron density inside the electronic surface. In Sec. IIID, we show the electronic stress tensordensityandtensiondensityoftheatoms. Finally, Sec.IVisdevotedtoourconclusion. II. THEORY AND CALCULATION METHODS In this paper, we analyze the electronic structure of atoms using quantities such as the electronic stress tensor density and the kinetic energy density. They are based on the RQED theory [2] and we briefly describe them in this section. For other studies of quantum systems with the stress tensor in somewhat different contexts and definitions, we refer Refs. [35–52]. We also refer Refs. [53–55] regarding other definitions of kinetic energy density. The most basic quantity in RQED is the electronic stress tensor density operator τˆΠkl(x) e which is defined as [1]. τˆΠkl(x) = i(cid:126)c (cid:20)ψˆ¯(x)γlDˆ (x)ψˆ(x)−(cid:16)Dˆ (x)ψˆ(x)(cid:17)†γ0γlψˆ(x)(cid:21), (1) e 2 ek ek ˆ where ψ(x) is the four-component Dirac field operator for electrons with the spacetime co- ordinate x = (ct,(cid:126)r). c denotes the speed of light in vacuum, (cid:126) the reduced Planck constant, 3 andγµ (µ =0-3)thegammamatrices. ThedaggerasasuperscriptisusedtoexpressHermite conjugate, and ψˆ¯(x) ≡ ψˆ†(x)γ0. The Latin letter indices like k and l express space coordi- ˆ nates from 1 to 3, and repeated indices implies a summation over 1 to 3. D (x) is the gauge ek covariant derivative and it is defined as Dˆ (x) = ∂ +iZeeAˆ (x), where Z = −1 and Aˆ (x) ek k (cid:126)c k e k ˆ (cid:126) is the vector potential of the photon field operator in the Coulomb gauge (divA(x) = 0). The important property of this quantity is that the time derivative of the electronic kinetic (cid:18) (cid:19) ˆ ˆ (cid:16)ˆ (cid:17)† momentum density operator Π(cid:126) (x) = 1 i(cid:126)ψˆ†(x)D(cid:126) (x)ψˆ(x)−i(cid:126) D(cid:126) (x)ψˆ(x) ·ψˆ(x) is e 2 e e ˆ (cid:126) expressed by the sum of the Lorentz force density operator L (x) and the tension density e operator (cid:126)τˆΠ(x), which is the divergence of τˆΠkl(x). Namely, e e ∂ ˆ ˆ Π(cid:126) (x) = L(cid:126) (x)+(cid:126)τˆΠ(x), (2) ∂t e e e ˆ ˆ 1ˆ ˆ L(cid:126) (x) = E(cid:126)(x)ρˆ (x)+ (cid:126)j (x)×B(cid:126)(x), (3) e e e c τˆΠk(x) = ∂ τˆΠkl(x) (4) e l e (cid:34) = i(cid:126)c (cid:16)Dˆ (x)ψˆ(x)(cid:17)†γ0γl ·Dˆ (x)ψˆ(x)+ψˆ¯(x)γlDˆ (x)Dˆ (x)ψˆ(x) el ek ek el 2 (cid:35) (cid:16) (cid:17)† (cid:16) (cid:17)† − Dˆ (x)Dˆ (x)ψˆ(x) γ0γl ·ψˆ(x)− Dˆ (x)ψˆ(x) γ0γl ·Dˆ (x)ψˆ(x) ek el ek el 1 (cid:16)ˆ (cid:17)k − (cid:126)j (x)×B(cid:126)(x) (5) e c ˆ where ρˆ (x) and(cid:126)j (x) are the electronic charge density operator and charge current density e e ˆ ˆ (cid:126) (cid:126) operator respectively, and E(x) and B(x) denote the electric field operator and magnetic field operator respectively. As we study nonrelativistic systems in this paper, we approximate the expressions above in the framework of the primary RQED approximation [4, 5]. In this approximation, the smallcomponentsofthefour-componentelectronfieldareexpressedbythelargecomponents as ψˆ (x) ≈ − 1 i(cid:126)σkD ψˆ (x), where m is the electron mass, and the spin-dependent terms S 2mc k L are ignored. Then, we take the expectation value of Eq. (2) with respect to the stationary state of the electrostatic Hamiltonian. This leads to the equilibrium equation as 0 = (cid:104)Lˆk(x)(cid:105)+(cid:104)τˆSk(x)(cid:105) = (cid:104)Lˆk(x)(cid:105)+∂ (cid:104)τˆSkl(x)(cid:105), (6) e e e l e 4 which shows the balance between electromagnetic force and quantum field force at each point in space. Since this expresses the fact that the latter force keeps the electrons in the stationary bound state in atomic and molecular systems, we can study them from the viewpoint of quantum field theory by using the stress tensor density and tension density. We express (cid:104)τˆSk(x)(cid:105) and (cid:104)τˆSkl(x)(cid:105) respectively τSk((cid:126)r) and τSkl((cid:126)r) for simplicity. Note that, e e e e as we consider stationary state, we write only spatial coordinate (cid:126)r. Writing explicitly, (cid:34) (cid:126)2 (cid:88) ∂2ψ ((cid:126)r) ∂ψ∗((cid:126)r)∂ψ ((cid:126)r) τSkl((cid:126)r) = ν ψ∗((cid:126)r) i − i i e 4m i i ∂xk∂xl ∂xk ∂xl i (cid:35) ∂2ψ∗((cid:126)r) ∂ψ∗((cid:126)r)∂ψ ((cid:126)r) + i ψ ((cid:126)r)− i i , (7) ∂xk∂xl i ∂xl ∂xk τSk((cid:126)r) = ∂ τSkl((cid:126)r) e l e (cid:34) (cid:126)2 (cid:88) ∂∆ψ ((cid:126)r) ∂ψ∗((cid:126)r) = ν ψ∗((cid:126)r) i − i ∆ψ ((cid:126)r) 4m i i ∂xk ∂xk i i (cid:35) ∂∆ψ∗((cid:126)r) ∂ψ ((cid:126)r) + i ψ ((cid:126)r)−∆ψ∗((cid:126)r) i , (8) ∂xk i i ∂xk where ψ ((cid:126)r) is the ith natural orbital and ν is its occupation number. ∆ denotes the i i Laplacian, ∆ ≡ (cid:80)3 (∂/∂xk)2. The eigenvalue of the symmetric tensor ↔τS is the principal k=1 e stress and the eigenvector is the principal axis. We denote the eigenvalues as τS11((cid:126)r) ≤ e τS22((cid:126)r) ≤ τS33((cid:126)r). e e We note that the tension density in the form of Eq. (8) is same as what is called quantum force density in Refs. [44, 56]. Then, the Ehrenfest force field used in Refs. [44, 56–58] (and the force density in Ref. [36], Eq. (24)) is same as the tension density with the minus sign in the stationary state. Another note is on the ambiguity of the stress tensor density. It is not defined uniquely since mathematically any tensor whose divergence is zero can be added to. Our stress tensor density Eq. (7) is same as the one in Ref. [36], Eq. (22). The difference in the definitions and approximations are discussed recently in Refs. [48–52]. We advocate the use of Eq. (7) since it comes from the stress tensor density operator Eq. (1), which is a minimal combination respecting the Lorentz covariance, gauge invariance and hermiticity. Moreover, this definition turns out to be phenomenologically useful as shown by our works mentioned in the previous section. 5 AnotherimportantquantityintheRQEDistheelectronickineticenergydensityoperator defined as [1], (cid:126)2 1 (cid:18) ˆ (cid:16)ˆ (cid:17)† (cid:19) Tˆ (x) = − · ψˆ†(x)D(cid:126)2(x)ψˆ(x)+ D(cid:126)2(x)ψˆ(x) ·ψˆ(x) . (9) e 2m 2 e e As is done for the electronic stress tensor density operator, we apply the primary RQED approximation to Eq. (9) and take the expectation value with respect to the stationary state of the electrostatic Hamiltonian. Then, we obtain the definition for the electronic kinetic energy density as (cid:126)2 (cid:88) n ((cid:126)r) = − ν [ψ∗((cid:126)r)∆ψ ((cid:126)r)+∆ψ∗((cid:126)r)·ψ ((cid:126)r)]. (10) Te 4m i i i i i i Notethatourdefinitionoftheelectronickineticenergydensityisnotpositive-definite. Using n ((cid:126)r), we can divide the whole space into three types of region as follows [1]: Te R = {(cid:126)r|n ((cid:126)r) > 0} : electronic drop region (11) D Te S = {(cid:126)r|n ((cid:126)r) = 0} : electronic interface (12) Te R = {(cid:126)r|n ((cid:126)r) < 0} : electronic atmosphere region (13) A Te In R , the electronic drop region, the classically allowed motion of electron is guaranteed D and the electron density is amply accumulated. In R , the electronic atmosphere region, the A motionofelectronisclassicallyforbiddenandtheelectrondensityisdriedup. Theboundary S between R and R is called the electronic interface and corresponds to a turning point. D A The outermost S can give a clear image of the intrinsic shape of atoms and molecules and is, therefore, an important region in particular. Thus, in our method, the outermost electronic interface defines the shape of atoms and molecules. Wenotethatitisfrequentlydefinedbytheisosurfaceoftheelectrondensity. The values of isosurface like 0.001 a.u. and 0.002 a.u. are proposed in Refs. [59–61], and there is some arbitrariness. In contrast, since our definition uses the zero isosurface of non-positive- definite quantity n ((cid:126)r), there is no such arbitrariness. It is true that there is ambiguity Te regarding the definition of the kinetic energy density [53–55], which is sometimes defined as a positive definite quantity. Our definition comes from the electronic kinetic energy density 6 operator, Eq. (9), whose definition is in turn motivated by the relativistic energy dispersion relation between the energy E and momentum p: E = (cid:112)(pc)2 +(mc2)2 ≈ mc2 + p2 . We 2m take the kinetic energy part p2 , replace p by i(cid:126)Dˆ as a usual quantization rule under 2m k ek the existence of the electromagnetic field, and construct a field operator by sandwiching between ψˆ†(x) and ψˆ(x). We then make the field operator Hermitian as Eq. (9) by adding the Hermitian conjugate and divide by two. We advocate the use of definition (10) since it comes from such a field theoretic construction. In the end of this section, we summarize our computational setups. The electronic struc- tures used in this paper are obtained by the Gaussian 09 [62]. The computation is per- formed by the unrestricted Hartree-Fock (UHF) method using the cc-pV5Z basis set [63]. To compute the aforementioned quantities such as Eqs. (7), (8) and (10) from the electronic structure data, we use the QEDynamics package [64] developed in our group. III. RESULTS AND DISCUSSION A. Shape of the electronic interface of atoms In Fig. 1, we show the electronic interface, S, for atoms of elements in the first, second and third periods. Spin multiplicities are chosen to give their ground states, whose term symbols are shown in Table I [65]. In the case of open-shell systems, the expectation values ofS2 arereportedinTableI,showingthatspincontaminationsaresmall. Asfortheelectron configuration regarding the p orbitals, we choose them to be occupied in the alphabetical order p , p , p . For instance, the configuration of B is 1s22s22p1, C is 1s22s22p12p1 and x y z x x y so on. In each panel, the atomic nucleus is placed at the origin. For all of these elements, the S’s have axial symmetry, and some of them have spherical symmetry. We omit the information of the value of the kinetic energy density in the figure, but we can recognize its sign by noting that the region in the neighborhood of the nucleus should have positive kinetic energy density, that is to say R . Going outward, R and R appear alternately D A D bounded by S and we should have R at infinity. In Fig. 2, we show the cross sections of A S’s to inspect their detailed structures. Again, in each panel, the atomic nucleus is placed at the origin. Since we plot to make the symmetry axis coincides with the horizontal axis, the axes may differ from panel to panel. For example, in the panel of B, the horizontal axis 7 corresponds to x-axis, but it corresponds to z-axis in the panel of C. We can see that the shapes found in Figs. 1 and 2 have the symmetry which is expected from the electron configurations of each atom. Namely, while the S’s of H, He, Li, Be, N, Ne, Na, Mg, P and Ar are spherically symmetric, those of the other elements, B, C, O, F, Al, Si, S and Cl are only axially symmetric. In detail, B, O, Al and S are axially symmetric with respect to the x-axis, reflecting the configuration of the outermost p subshell: np1 for x B and Al, and np2np1np1 for O and S. As for C, F, Si and Cl, they are axially symmetric x y z with respect to the z-axis due to the configuration np1np1 for C and Si, and np2np2np1 for F x y x y z and Cl. Also, the periodicity regarding the size of S is manifestly shown: the size increases as the period increase in the same row, and decreases as the atomic number increases in the same period. The common feature for all the elements is that they clearly have outermost S which is homeomorphic to 2-sphere (two-dimensional surface of a ball). We denote the outermost S as S . We next notice that, while H, He and Ne have only one S (only S ), the others outer outer have multiple S’s. As for those with multiple S’s, most of them have two S’s inside S , outer and in addition, they are both homeomorphic to 2-sphere and the larger S encloses the smaller ones. When such a structure is found, we denote the innermost S as S . In every inner case, since these two S’s are very close to each other, it seems that a relatively thin spherical shell region of R is formed inside S . In other words, the region inside S appears A outer outer to be divided into two R ’s by the nearly spherical-shell-shaped R . Exceptions are O and D A F. As for O, it has a single S which is homeomorphic to 2-torus (a solid-torus-shaped R ) A inside S . As for F, it has two S’s inside S but they do not form a shell-like region outer outer and just two disconnected R ’s are formed. (Actually, these features are so fine that it is A very hard to see in Fig 1. They can be more easily imagined from Fig. 2 by rotating the cross section with respect to horizontal axis.) Thus, for both cases, there is only one connected region of R , and S is not formed. We check whether these features for O and F may D inner change if we use wavefunctions computed by the restricted open-shell HF (ROHF) method. We find that the ROHF gives visually undistinguishable S’s for O and F. We discuss why they form such R ’s below. A On looking at these structures, one might associate the pattern of S as the electron shell structure in atoms. It is tempting to associate R which is divided by thin R as the shell D A structure. The single R seen in H and He can be associated with K shell. From Li to D 8 N, they have two R ’s which can be associated with K and L shells. However, this breaks D down at O, F and Ne, which have only one R as pointed out above. Moreover, from Na to D Ar, which should have K, L and M shells, they only have two R ’s. Therefore, S is not a D good indicator of the atomic shell structure. Instead, it may be better to consider two R ’s D separated by the thin R region as the manifestation of “core” and “valence” regions in A some sense. In particular, the core shell can be defined as R enclosed by S if it exists. D inner The irregularities concerning O, F and Ne are interpreted as follows. As for O, if there were only three electrons in the 2p subshell as 2p12p12p1, the R would be shell-shaped region x y z A like the one in N. Since there exists one more 2p electron, the probability of the electron x in the core region is higher so that shell-shaped R is deformed to be solid-torus shape as if A penetrated by the 2p orbital. Note that the torus’s hole is opened in the direction of x-axis, x which is consistent with the direction of 2p orbital. As for F, one more electron in the 2p x y orbital is added and it again penetrates in the core region. Then, the torus-shaped R is A split in the y-direction, and two R ’s remain on the z-axis. Finally, in Ne, the addition of A one more electron in the 2p orbital completely erases these remaining R ’s. z A B. Size of the electronic interface and comparison with atomic and ionic radii In this section, we quantify the size of S and S in the following manner. When outer inner S is spherical, the size can be readily defined by its radius. However, as we have seen in Sec.IIIA,someoftheelementshavenon-sphericalS. Inthatcase, wereportthelongestand shortest distances between nucleus and S, respectively as the upper and lower limits. Also, as the central value, we report an averaged distance between nucleus and S by computing the radius of the ball whose volume is equal to the S’s volume. As for O and F, for which S is not defined, we just report the shortest distance between nucleus and R inside inner A S . The results are summarized in Table I and Fig. 3. We can confirm the periodicity outer regarding the size of S and S in the figure. We note that, for S , the deviation outer inner inner from sphere is so small that we cannot recognize the error bars. As for O and F, the same quantities are computed using the ROHF wavefunctions. We find that S and S of outer inner O are respectively 0.9610+0.0446 and [0.1533], and those of F are respectively 0.8623+0.0159 −0.0313 −0.0451 and [0.1378], which are almost identical to the UHF results in Table I. At this stage, it is instructive to examine S of the cations which are constructed by 9 removing all the electrons in the outermost shell, for example, B3+, C4+, N5+ and so on. As a result, we find that all these cations have only one S (only S ) which is spherical, outer as expected. The S ’s radii of the cations are summarized in Table I and Fig. 3. Then, outer we see that S of each atom is almost same as the S of the corresponding cation. inner outer This is another reason that we wish to interpret R bounded by S as a core region. D inner (In addition, we have observed in our past works that such core regions defined by S are preserved in molecules. For example, Figs. 2, 3, 7 and 8 of Ref. [32] for Li and Fig. 2 of Ref. [33] for C.) As for O and F, although they do not have S , the sizes of R regions inner A inside their S are very close to S of O+6 and F+7. This supports our interpretation outer outer for the irregular shapes of R regions in O and F, that they are caused by the penetration A of electrons in 2p orbitals into their core regions, as argued in the end of Sec. IIIA. Next, in Fig. 3, we plot atomic radii from Ref. [66] and ionic radii from Ref. [67], along with the size of S and S as determined above. We quote two types of ionic radii, outer inner crystal radii and univalent radii, which appear in Ref. [67]. The crystal radii of multivalent ions are defined so that the sum of two radii is equal to the actual equilibrium interionic distance in a crystal. For example, the univalent radius is applied to Mg2+F− while the 2 crystal radius to Mg2+O2−. We find that the atomic radii are comparable to the size of S , in a sense that they have similar periodicity and differ by a factor of two at most. outer Similarly,wecansaythattheionicradiiarecomparabletothesizeofS . However,taking inner a closer look, the atomic radii are always smaller than the size of S . The univalent ionic outer radii are always larger than the size of S . Although some of the crystal ionic radii inner are very close to the size of S , since the slopes with respect to the nuclear charge are inner different, the agreement should be considered as a mere coincidence. Such a disagreement between the atomic radii and S ’s size, and the one between outer the ionic radii and S ’s size, are not surprising. This is because, while the atomic/ionic inner radii are defined to give equilibrium internuclear distance in compounds and crystals on summing two radii, the sizes of S and S are defined for isolated atoms. As for outer inner the atomic radii, we can consider that the formation of covalent bond just begins when the S ’s of two atoms touch each other (this moment is called the “electronic transition state” outer [1]). Then, after the bond formation, at the equilibrium, the interatomic distance is shorter than the sum of the radii of S of the two atoms. Therefore, the formation of chemical outer bonding makes the atomic radii smaller than the S ’s sizes. As for the ionic radii, let outer 10