Topological Analysis of the Electron Density Illustrating the Stabilizing Interactions in Some Basic Solid-state Systems

The solid-state structures of Sodium (Na), Titanium (Ti), Diamond and Graphite, Sodium Chloride (NaCl), Magnesium Oxide (MgO), Cadmium (II) Iodide (CdI2) and Zirconium Chloride (ZrCl) have been explored in details using computational electron density methods; the full-potential linearized augmented plane wave (FPLAPW) method plus local orbital (FPLAPW+lo) embodied in the WIEN2k package code. Topological analysis of their DFT-computed electron densities in tandem with Bader’s Atoms in Molecules (AIM) theory reveals a plethora of stabilizing interactions some of which are really strong. Na and Ti metals reveal only metallic bonding, diamond and graphite show covalent bonding between the carbon-atoms. In addition, there exist Van der Waals forces between the carbon-atoms on adjacent planes in the graphene sheets. NaCl and MgO exhibit electrostatic interactions between the metals (Na, Mg) and non-metals (Cl, O) respectively. Furthermore, there exist Van der Waals interactions between Cl and O atoms. CdI2 and ZrCl both show ionic and Van der Waals forces between the atoms. ZrCl exhibit metallic bonding and NNMs between the Zr-atoms, which are absent in CdI2 due to longer Cd-Cd bond distances. Analyses of the electron density flatness (f), charge transfer index (c) and molecularity (μ) were computed. It is observed that the different types of interactions increase with complexity of the solid-state structures. Finally, non-nuclear maxima (NNM) were identified for the first time in heteroatomic solid-state systems.


Introduction
Atomic packing is of paramount importance in classifying the different types of interactions (bonding) within crystals (complexes) [1]. These crystals may be finite or extend indefinitely in one, two, or three dimensions and stabilized by ionic, covalent, metallic, Van der Waals or hydrogen bonds. Recognition and identification of the stabilizing interactions provides the basis for a broad geometrical classification of crystal structures. However, approximations of bonds in pure-type crystals (ionic or covalent) are rare, and moreover in most crystals there exist bonds of different types. In some cases, numerous intermediates have to be recognized, and this makes classification of crystals based on bond types to be complicated and incomprehensive. It is therefore necessary investigate, analyze and discuss the nature and types of bonds in crystals without having prejudged the issue by classifying them as ionic, for example [1]. Furthermore, understanding the electronic and atomistic principles of solids is essential for practical use. The ways in which atoms are arranged and interact with one another determines the fundamental characteristics of the solid state systems (crystals) and their possible changes tell us how we could modify the materials for various applications. In this paper, a systematic investigation on the bond types that stabilize some solid-state systems is presented using Bader's Atoms in Molecules (AIMs) theory. The chemical and physical characteristics of molecules in general and solids in particular, are based on the nature of the interactions that exist between the atoms. The theory of 'atoms in molecules' AIM [2,3] can be seen as an ongoing research that aims at extracting and interpreting chemical information from modern quantum mechanical techniques (ab initio wave function) in solids. It was earmarked or propounded in the early 1970s by Richard Bader and co-worker [2], and is currently being used by more than 70 laboratories worldwide in areas of surface science, organometallic chemistry, life science, solid-state physics and chemistry, drug design, physical-organic chemistry, crystallography, etc. The theory provides a simple, rigorous and elegant way of analyzing and interpreting the various types of stabilizing interactions in solid materials (maintaining atoms individuality) from the electron density and its laplacian [2,[4][5][6][7].
The purpose of this paper is to investigate and interpret the types of stabilizing interactions in some solid-state systems ranging from simple metals (Na, 1, Ti, 2), nonmetals (diamond, 3 and graphite, 4), simple binary ionic compounds (NaCl, 5, MgO, 6), and complex binary systems (CdI 2 , 7, ZrCl, 8). While metals and non-metals are expected to exhibit a single type of bonding, binary systems show more than one type of interaction. This is as a result of the complexity in their crystal structures. We observe that as the complexity of the crystal structure increases, they will require more than one type of bonding to effectively stabilize their structure. As a result, simple metals and non-metals are expected to exhibit only metallic and covalent bonding respectively. The more complex ionic crystals 5 and 6, are expected to exhibit more than one type of interaction; a strong closed-shell ionic interaction (electrostatic bonding), which result in a complete transfer of electrons from the metal to the nonmetal and weak van der Waals interactions between the non-metallic moieties that depend on their distances. In addition, we would expect the existence of non-nuclear maxima (NNM) of electron density in some of these systems. NNMs are local maxima of electron density that appears between two or more nuclei. They behave as pseudoatoms, trapping most of the valence electrons of the metallic atoms and occur at interstitial positions other than the atomic nuclei and are topologically indistinguishable from the nuclear maxima. This phenomenon is common and has been observed in homopolar diatomics [8][9][10][11][12][13][14]. No previous work on the existence of these NNM's in heteroatomic systems have been reported in literature [12]. Finally, the types of interactions in these systems will also be analyzed and interpreted based on the three indexes vis; electron density flatness, charge transfer index and molecularity [9].

Methodology
We used the full-potential linearized augmented plane wave (FPLAPW) method plus local orbital (FPLAPW+lo) embodied in the WIEN2k package code [15]. No shape approximations are applied for charge density or potential. The exchange-correlation effects are treated in the density functional theory (DFT) within the FPLAPW formalism, using the generalized gradient approximation (GGA) together with Perdew and Wang functional (PBE96) [16]. The crystal unit cell used in this approach is partitioned into non-overlap atomic spheres (muffin tins) and interstitial region. The radii of the muffin tin spheres are constrained by the requirement that they are non-overlapping and that the core states do not significantly spill into the interstitial region. Inside each muffin tin, orbitals are described as radial functions times spherical harmonics and a fully relativistic treatment is used, based on the work by Desclaux [17]. The interstitial region is described using plane waves and a scalar relativistic treatment is applied [18]. The wave functions of both regions are set to match in value and slope at the muffin tin boundaries, but a perfect match would require the inclusion of spherical harmonics to infinite order within each atomic sphere and a truncated series is instead used. Care must be taken to ensure that any remaining discontinuity in the density or its slope is not seen by the topological algorithms. Keeping the muffin tin radii smaller than the shortest distance from a nucleus to any critical point of the electron density is recommended. The following muffin tin sphere radii were used for the calculations: Na (2.5 au) for 1, Ti (2.5 au) for 2, C (1.5 au) for 3 and C (1.3 au) for 4, Na (2.5 au) and Cl (2.5 au) for 5, Mg (1.88 au) and O (1.88 au) for 6, Cd (2.5 au) and I (2.5 au) for 7, Zr (1.5 au) and Cl (1.8 au) for 8. The cut-off parameter R mt ×K max for limiting the number of the plane waves is equal to 7, where R mt is the smallest of all atomic sphere radii and K max is the largest reciprocal lattice vector used in plane wave expansion. For the Brillouin zone integration, we used K-points mesh values of 10×10×10 in 1, 10×10×10 in 2, 14×14×14 in 3, 10×10×10 in 4, 11×11×11 in 5, 9×9×9 in 6, 16×16×25 in 7, and 34×34×35 in 8 for the first Brillouin zone. The self-consistency was achieved when the total energy was found to be stable within 10 -4 Ry. The topological analysis of the electron density of these solidstate systems were done using CRITIC [19], a quantum theory of atoms in molecules (AIM) [2] topological code.

Analysis of Sodium (1)
Sodium crystallises in a cubic (bcc) unit cell with space group (229) Im3m and cell parameters a = b = c = 4.2906 Å, α = β = γ = 90. The topological parameters after full optimization of the structure are shown in Table 1. We observe metallic bonding between the Na-atoms in 1, which is represented by the very low values of the electron density, ρ b (r) (0.031 e/Å 3 ) and its positive laplacian, ∇ 2 ρ b (r) (0.009 e/Å 5 ) approaching zero. Electrons are seen to be localized on the atoms as shown on the laplacian plot, Figure 1.

Analysis of Titanium (2)
Titanium exist in the cubic space group F m -3m (225) with unit cell parameters a = b =c = 4.06 Å, α = β = γ = 90 [20]. Analysis of 2, revealed significant differences in interaction when compared to 1. There exists stronger metallic bonding between the Ti-atoms in 2 when compared to Na-atoms in 1. This indicates that Ti is more metallic that Na with observed polarization on the Ti-atoms. These results are further confirmed by the higher values of the electron density, ρ b (r) (0.201 e/Å 3 ) and its positive laplacian, ∇ 2 ρ b (r) (0.080 e/Å 5 ). See results in Table 1.

Analysis of Diamond (3) and Graphite (4)
Diamond crystallizes in the cubic lattice with space group F d -3s (227). The unit cell parameters are a = b = c = 3.56672 Å, α = β = γ = 90 [21]. Graphite on the other hand crystallizes in the hexagonal lattice with P 63/mmc (194) space group and unit cell parameters a = b = 2.464 Å, c = 6.711 Å [22]. In diamond, each C-atom is linked to four equidistant (1.5 Å) neighbours and extends throughout the whole crystal while in graphite, the layers are separated by a distance of 3.4 Å, with C-C bond lengths of 1.4 Å. This are in agreement with literature values [1]. Results from the topological parameters reveal that, the atoms in 3 are held by strong C-C covalent bonds. This is confirmed by the high electron density, ρ b (r) (1.623 e/Å 3 ) and negative laplacian, ∇ 2 ρ b (r) (-16.53 e/Å 5 ) values. Conversely, 4 possesses two types of bonding viz; covalent bonding between the C-atoms in the same plane and weak Van der Waals interaction between the graphene sheets (C-atoms of different layers). These results are also confirmed by the electron density and laplacian values. The C-C covalent bonds have high electron density, ρ b (r) (2.004 e/Å 3 ) and negative laplacian, ∇ 2 ρ b (r) (-25.2 e/Å 5 ) values. The van der Waals interactions between the graphene layers are shown by the low values of the electron density, ρ b (r) (0.037 e/Å 3 ) and positive laplacian, ∇ 2 ρ b (r) (0.376 e/Å 5 ). Changes in the interactions can also be accounted for by the C-C distances. Furthermore, the ellipticity, ε, indicates the degree of deviation from sigma symmetry. This is confirmed by the spherically symmetric C-C single bonds in 3 (ε = 0.00) and the non-symmetric delocalisation in 4 (ε = 0.16), with ellipticity lying between that of a purely single bond (ε = 0.00) and double bond (ε = 0.34). The topological parameters are shown in Table 1.

Analysis of Sodium Chloride (5)
Binary ionic crystals such as 5 and 6 crystallizes in the cubic lattice with F m-3m (225) space group. The unit cell parameters for 5 are, a = b = c = 5.6400 Å, α = β = γ = 90 [23] and for 6, a = b = c = 4.217 Å [24]. These structures show an alternate arrangements of the atoms with regular octahedral coordination of both atoms. The crystal structure of NaCl is shown in figure 3. The results in Table 1 indicate ionic or electrovalent bonding between Na + and Clions in 5, Mg 2+ and O 2ions in 6 and Van der Waals interactions between the Cl-atoms and O-atoms in 5 and 6 respectively as expected. The ionic interactions between atoms in these systems is illustrated by the complete transfer of electrons from Na to Cl and from Mg to O in 5 and 6 respectively ( Table 1) (Table 1). Figure 3 represents the structure and laplacian plot in a plane of NaCl.

Analysis of Cadmium Iodide (7)
Cadmium iodide crystallizes in the trigonal/rhombohedral lattice with space group P -3 m1 (164) and with unit cell parameters a = b = 4.244 Å, α = β = 90, γ = 120 [25]. The structure of 7 is more complex compared to 5 and 6 and composed of alternate layers of Cd and I 2 . The complexity of the structure will therefore, imply that they show different types of interactions. The crystal structure of 7 is shown in Figure 4. Our topological analysis revealed electrostatic bonding between Cd 2+ and Iions and london dispersion forces between the I-atoms. This is confirmed by the transfer of electrons from Cd to I 2 , Exhibiting low values of electron density, ρ b (r) (0.220 e/Å 3 ) and positive laplacian, ∇ 2 ρ b (r) (0.747 e/Å 5 ). On the other hand, the weak interactions between the I-atoms is shown by the low values of the electron density, ρ b (r) (0.031 e/Å 3 ) and positive laplacian, ∇ 2 ρ b (r) (0.245 e/Å 5 ). Surprisingly, no metallic bonding was observed between the Cd-atoms. This is probably due to the long Cd-Cd bond distances (> 3.01 Å) between the Cd-atoms on adjacent layers.

Analysis of Zirconium Chloride (8)
The complexity of our systems were further extended to structure 8, exhibiting different types of interactions. It is composed of alternating layers of Zr and Cl. Structure 8 crystallizes in the R -3m H (166) space group with unit cell coordinates a = b =3.424 Å, α = β = 90, γ = 120 [26]. It exhibits three (3) different types of interactions viz; ionic bonding between Zr + and Clions, london dispersion forces between the Clatoms and metallic bonding between the Zratoms. These interactions are confirmed by Their topological parameters shown in Table 1. Ionic bonding is observed between Zr + and Clions exhibiting electron density, ρ b (r) (0.279 e/Å 3 ) and positive laplacian, ∇ 2 ρ b (r) (2.830 e/Å 5 ) values. London dispersion forces are seen between the Cl-atoms on adjacent layers. This is represented by the low values of the electron density, ρ b (r) (0.010 e/Å 3 ) and positive laplacian, ∇ 2 ρ b (r) (0.114 e/Å 5 ). In addition, there exists metallic bonding between the Zr-atoms in the same plane. This is confirmed by the electron density, ρ b (r) (0.196 e/Å 3 ) and laplacian, ∇ 2 ρ b (r) (0.138 e/Å 5 ) values. It should be noted that during the formation of electrostatic bonding, the electrons are not completely transferred from Zr to Cl. Therefore, charge of 0.134 e/Å 3 is concentrated between the Zr-atoms as shown in Figure 5. This was identified as nonnuclear maxima (NNM) or non-nuclear attractors (NNAs) between the Zr-atoms. The phenomenon was first discovered by Besnainou in 1955 [27], and later in the early 1980s by Richard Bader [28]. This was further supported by several reports in the 1980s and early 1990s [29][30][31][32][33][34][35][36]. A final consensus was arrived at in 1999 that NNMs are stable features on the theoretical electron density in some systems [12,13]. NNMs are local maxima of electron density that appears between two or more nuclei, and behave as pseudoatoms, trapping most of the valence electrons of the metallic atoms. They occur at interstitial positions other than the atomic nuclei and are topologically indistinguishable from the nuclear maxima with no interactions observed between them. This has been reported in most metals but not in heteroatomics [12]. Therefore, the existence of NNMs in heteroatomics (eg ZrCl) provides a unique feature in these types of solid-state systems.

Classification of the Interactions Based on the Three Indexes
Previous report from Martin Pendas (2002) classified covalent, ionic and metallic bonding in simple solids using the electron density. Three parameters (flatness, charge transfer index and molecularity) were used in the classification scheme [7]. The results are presented in Table  2.

Conclusion
The topologies of the electron density for the above solids (1-8) have been fully analyzed and interpreted. It is seen that the types of stabilizing interactions in these systems increase with increase in the complexity of the solid-state structures. Sodium and Titanium metals exhibit only metallic bonding as expected, while diamond and graphite show covalent interaction between the C-atoms. In addition, Van der Waals interactions were observed between the C-atoms on adjacent planes in graphene sheets. Binary ionic crystals such as sodium chloride, 5 and magnesium oxide, 6 exhibit electrostatic interactions between the metals (Na, Mg) and the non-metal (Cl, O) respectively as expected. We could also identify some weak van der Waals interaction respectively between the Cl and O-atoms due to the complexity of the NaCl and MgO structures. Furthermore, Cadmium iodide (CdI 2 ) exhibits ionic and Van der Waals interactions between the atoms. Surprisingly, no metallic bonding was observed between the Cd-atoms due the long Cd-Cd distances. Zirconium chloride (ZrCl) exhibits ionic, van der Waals, and metallic bonding between the atoms (Table 1). In Addition, there exist NNMs within the ZrCl structure with no interactions observed between the NNMs and Zr or NNMs and Cl. This is the first time NNMs (NNAs) are reported in heteroatomic systems. Finally, our systems were analyzed based on the three indexes proposed by Pendas and the results obtained were agreeable within less 5% error margin.

Author Contributions
J. T conceived, put together, and performed the DFT calculations described in this manuscript with assistance from F. C. The manuscript was written with contributions from all authors.

Funding Agencies
We are thankful to the Natural Sciences and Engineering Research Council (NSERC) of Canada and its H2CAN network for financial support.