Topology vs. thermodynamics in chemical reactions: the instability of PH₅

Abstract

The topological approach, based on Bader theory, is compared to the common thermodynamical methodology to study chemical reactivity. It is shown how the former indeed has numerous advantages and provides a more detailed description with respect to the latter about the course of the reaction. The comparison between the two approaches is performed by considering a classical reaction, i.e. the decomposition of PX₅ (X = H, F). The topological investigation was supported by using different state-of-the-art topological tools, such as the source function, Espinosa indexes, delocalisation indexes, and domain-averaged Fermi hole analysis. Furthermore, in this work a new topological descriptor, the Bader energy density, P_BADER, is introduced and applied to the study case. For the first time since Bader theory was introduced, the distribution of atomic energies in the atomic basins was analysed in detail and used to explain the chemical reactivity a priori.

---

Introduction

The last decades have seen the development of modern theories to study charge density.¹ Various theories, including the quantum theory of atoms in molecules (QTAIM),² have found important applications in the understanding of chemical, physical and biological phenomena. The description of existing correlations between electronic structure and chemical reactivity is a research field of particular interest and relevance.^3–9 In this field the main objective is to predict and to describe the reaction path of an arbitrary process a priori, simply through analysis of the charge density. In some reactions involving “heavy” elements (e.g. P, S, etc.), they may change their valence shell violating the so-called “octet rule:”¹⁰ an expansion of the atom's valence shell is the cause of filling empty d-orbitals with energy close to those of valence orbitals. This description – introduced by Pauling^3,11 under the name of hypervalency – is one of the most considered models employed to illustrate and interpret the electronic structure of “heavy” elements. Through studies based on the QTAIM approach, the theory of “hypervalence” was proven to incorrectly describe the electronic structure of different compounds, such as H₂SO₄, PF₅ and sulphonamide groups within oxicams.^12–14 In these studies it was shown that there is no expansion of the valence shell with violation of the octet rule as previously proposed by Pauling.^3,11 Nevertheless the “hypervalency” theory is still widely evoked in many different chemistry text books, being the easiest way to give an apparently convincing description of chemical bonding, which does not require any knowledge about the quantum mechanics and electronic structure, as for Bader's theory. For example in 2015,¹⁵ Durrant took into account different Pauling- and Musher¹⁶-like models to examine the “hypervalency” of “heavy” elements in some compounds. These findings were compared with theoretical and experimental atomic charge maps in the framework of Bader's theory showing no fundamental differences in chemical bonding for hypervalent and non-hypervalent species. Nevertheless, hypervalency is associated with chemical instability, as well as a high degree of covalent rather than ionic bonding.

In this work we consider as a working example a model reaction involving a “heavy” atom that shows valence shell modification passing from “hypervalent” to non-“hypervalent” during the reaction; in particular we shall show our investigation on the well-known model decomposition reaction:

PX₅ = PX₃ + X₂, (1)

where X = H or a halogen. The quantum chemical rationalisation of phosphorus “hypervalency” still remains obscure. Phosphorus pentacoordinate compounds PX₅ are known and characterised using ligands, such as fluorine or chlorine atoms (X = F, Cl), whereas PH₅ is experimentally undiscovered due to its instability. The instability of PH₅ is supported by negative value of free energy change, ΔG^e_r < 0, for the decomposition reaction.¹⁷ Over the years there have been attempts to provide a theoretical rationalisation of such instability, and different explanations were proposed, often in contradiction^18–23 with each other. It seems that at a certain point the problem was abandoned, leaving the question unanswered, probably due to a lack of sufficiently powerful computational resources that could allow advanced calculations to be performed. We want to raise the problem again and attempt to provide a sound answer to it; this investigation is still worthy of interest and may shed light on the current topic of chemical reactivity in the presence of “hypervalent” atoms.

The Pauling model^3,11 describes phosphorus within PX₅ expanding its valence shell by filling its d-orbitals through dⁿsp³ hybridisation. In the case of PH₅ the superposition between s-orbital (belonging to H) and d-orbital (belonging to P) is avoided due to their high energy difference. This model was widely accepted until 1999, when Ponec et al.¹² showed (for a pentacoordinate compound of phosphorus, i.e. PF₅) that the inclusion of d-orbitals in the wavefunction is an artefact, being a result of polarisation functions added to the basis set (i.e. 6-311++G**). By means of domain-averaged Fermi hole (DAFH)¹² analysis, inside PF₅ (D_3h) the presence of a 3-centre-4-electron (3c-4e) bond between two axial fluorine atoms and the phosphorus atom, P–F_ax, was shown. The same theoretical framework may also be evoked to interpret the P–H_ax bond inside PH₅, which could be more likely interpreted as an interaction between a non-hybridised p-orbital of phosphorus and two s-orbitals of the apical hydrogen atoms. Therefore, hypothesising a 3c-4e bond within PH₅, the orbital superposition between s (of hydrogen) and p-orbitals (of phosphorous) might not be energetically avoided. However, Ponec et al. showed²⁴ that when the flexibility of the basis set is increased or electronic correlation is included, the 3c-4e bond in PF₅ should be interpreted as two very highly polarised 2-centre-2electron (2c-2e) bonds. At the moment, there is not a valid quantum-mechanical approach to clarify this question.

In the present work, we opt for a different approach. Instead of relying on the “classical” thermodynamic methodology, we attempt to rationalise the reactivity of PX₅ focusing our attention on the difference in topological features of the molecular systems involved in the decomposition reaction (1). In particular, our investigation develops within the framework of QTAIM,² by evaluating the source function (SF),^25–35 Espinosa indexes,³⁶ delocalisation indexes,³⁷δ(I,J), DAFH analysis,^12,24 and Bader energy density, P_BADER. The latter consists of a new topological descriptor here introduced and tested to explain the reactivity of PX₅ (vide postea). We finally show how a topological approach may have advantages and provide a more detailed picture with respect to a thermodynamic study in the investigation of chemical reactivity.

Theoretical framework

The introduction of a new topological tool was born from the idea that the total energy of an isolated system is conserved, even though during a reaction there is a partial energy transfer from reagents to products. At the same time also the total charge density for an isolated system is conserved and only the atomic charge distribution is modified upon partial charge density transfer from reagents to products. As clarified by Bader, “[⋯] its contribution to the total energy will be different and one may relate its change in energy to this change in its distribution of charge as induced by the new environment.”² Furthermore, the atomic basin has its own volume, wherein the charge density is enclosed, which is non-spherical due to the bond polarity and due to the conditions imposed by the Bader theory.^2,38 Thus, during a reaction, to take into account energy and charge transfer associated with each atomic basin, it is necessary for the latter to consider the variation of volume of each atomic basin. On the basis of these elementary considerations, the new topological descriptor P_BADER was introduced, consisting of the atomic basin electronic energy, E(Ω), per unit basin volume, V(Ω):

P_BADER(Ω) = E(Ω)/V(Ω) (2)

P_BADER may be, therefore, interpreted as an attempt to weigh the energy associated with the atomic basin accounting for its volume variation from one stationary point (e.g., reactant) to another (e.g., product).

The atomic basins are delimitated by a surface, S, the so called zero flux surface (ZFS),² defined by the geometrical locus of the points of the space, r̄_S, for which the scalar product of the gradient of charge density, ∇ρ, and the normal vector, n̂, is zero:

∇ρ(r̄)·n̂(r̄) = 0, ∀r̄ ∈ S(r̄_S) (3)

The ZFS is numerically approximated not only between two atomic basins, but also in the bound against a vacuum as the region wherein the charge density, ρ, is less than 10⁻³ e Bohr⁻³ considering that ρ tends to 0 when its distance from the nucleus tends to infinity. Under this approximation, Bader showed² that the accuracy of numerical integration is guaranteed by L(Ω) defined as:If L(Ω) falls in the range of 10⁻⁴ to 10⁻⁵, the error in the estimation of kinetic energy was demonstrated to be less than 0.4 kJ mol⁻¹. In the definition of P_BADER the same numerical approximation is assumed.

Within the framework of Bader theory,² the partition of molecular electron energy in atomic energies is carried out in observance of virial theorem; for an equilibrium geometry, the ratio γ = 〈V̂〉/〈T̂〉 = −2, is not satisfied. Thus, under these conditions, the expectation values of 〈T̂〉 and 〈V̂〉 must be multiplied by a factor (1 + γ) and (1 + 1/γ), respectively, in order to obtain the correct results E = −〈T̂〉 = 〈V̂〉/2.

The P_BADER is a simply intuitive tool, which allows the description of a reaction course.

Computational details

The molecular structures of PX₅, PX₃, and X₂ (X = H, F) were fully optimised in the gas phase at the following levels of theory: Hartree–Fock (HF), Møller–Plesset second order perturbation theory (MP2), coupled cluster (CC) for single and double excitation (CCSD), and quadratic configuration integration for single and double excitation (QCISD). These methods were coupled with the Dunning triple-ζ local basis set aug-cc-pVTZ. The levels of theory were chosen to account for the increase in electron correlation, which in the presence of F₂ was shown to play an important role.³⁹ The molecular symmetry of PX₅ molecules was the same, i.e. the symmetry point group D_3h with both elongated and non-elongated axial bond distances.¹² Associated with each molecular structure the PROAIM wave-functions⁴⁰ were calculated at each different level of theory for the subsequent Bader's analysis.

The calculation of delocalisation indexes, δ(I,J),³⁷ is based on the Fermi hole. Since in HF theory the Fermi hole was shown to be the sole source of correlation between the electrons,³⁷ the delocalisation indexes were computed using the PROAIM wave-functions obtained at the HF/aug-cc-pVTZ level.

For the same reason, the domain averaged Fermi hole (DAFH) analysis was also performed only for PH₅ and PF₅ molecular structures previously optimised at the HF/aug-cc-pVTZ level. The PROAIM wave functions⁴⁰ were, subsequently, computed using different basis sets, viz. 6-31G**, 6-311G**, and cc-pVTZ, without diffusion functions.

The full geometry optimisation in the gas phase with the associated PROAIM wave-functions was performed using Gaussian G09.D01.⁴¹ The Bader analysis was performed using a modified version of PROAIMV.^42,43 The PROAIM wave-functions for DAFH analysis were computed using Gaussian G03.C02.⁴⁴ The DAFH analysis was performed using the latest version of WinBader.⁴⁵

P _BADER is not implemented in any released program yet, but it is possible to calculate through a PROAIMV^42,43 calculation.

Results and discussion

Pentacoordinate phosphorus, PX₅, displays elongated (P–X_eq ≠ P–X_ax) and non-elongated (P–X_eq = P–X_ax) minima,¹² and they both were investigated with the same levels of theory, by looking at the BCP density and the source function percentage values (SF%). The latter reads:

SF%(r̄,Ω) = SF(r̄,Ω)/ρ(r̄)·100 (5)

The results gave the same qualitative picture for the two aforementioned minima (see Tables S3–S6, ESI†). Hereinafter we, therefore, present and discuss only the results for the elongated forms of PX₅, which is the form seen in nature, along with PX₃ and X₂.

Examining BCPs within PH₅ and PF₅ one finds that the charge density in axial BCPs is smaller than in equatorial BCPs with a difference in charge density of ca. 0.01 e Bohr⁻³ for CC and MP2, and ca. 0.11 e Bohr⁻³ for HF (Tables S3 and S4, ESI†). Nevertheless, the different charge density evaluated in BCPs is not particularly informative by itself. At the highest level of theory (viz., QCISD) the value of the Laplacian of charge density, ∇²ρ, calculated at BCPs along with other topological parameters according to Espinosa classes³⁶ (see Table 1), are also evaluated to examine the type of bonding interaction for all investigated systems. The topological parameters evaluated at BCPs of the PX₅ molecule showed different pictures for X = F with respect to X = H, and different numerical values for axial and equatorial BCPs: PH₅ shows a quite weak shared shell interaction,³⁶ as evident from the ratio between the absolute value of virial (V_b) and kinetic (G_b) energy at the BCP, Table 1, associated with a covalent bond.³⁶ Thus, the associated bond degree parameters, H_b/ρ_b were subsequently calculated to evaluate the interaction energy (see Table 1); PF₅ shows a transition zone,³⁶ which is commonly found in the case of metal bonds²⁶ or hydrogen bonds (H-bonds).³⁶

Table 1 Topological parameters at the BCP for the system involved in the dissociation reactions, investigated at the QCISD/Aug-cc-pVTZ level

System	Bond	ρ b (e Bohr^−3)	∇^2ρb^b (e Bohr^−5)	|Vb|^c (a.u.)	G b (a.u.)	H b (a.u.)	|Vb|/Gb	H b/ρb^f (a.u.)/(e Bohr^−3)
a Charge density. b Laplacian of charge density. c Virial energy. d Kinetic energy. e Hamiltonian. f Bond degree parameter.
PF5 (elongated)	P–Feq	0.19	1.15	0.56	0.42	−0.14	1.32	—
	P–Fax	0.18	0.94	0.50	0.37	−0.13	1.36	—
PH5 (elongated)	P–Heq	0.17	−0.06	0.33	0.16	−0.17	2.10	−1.04
	P–Hax	0.16	−0.11	0.30	0.14	−0.16	2.21	−1.05
PF3	P–F	0.17	0.98	0.48	0.36	−0.12	1.33	—
PH3	P–H	0.16	−0.01	0.32	0.16	−0.16	2.01	−1.00
F2	F–F	0.30	0.48	0.51	0.32	−0.20	1.62	—
H2	H–H	0.27	−1.19	0.32	0.01	−0.31	29.85	−1.15

---

Applying Espinosa's classification criteria,^9,36 PF₃ and F₂ contain interactions classified as transition zones, whereas PH₃ and H₂ contain interactions classified as shared shell. The difference between PF₅ and PH₅ is furthermore confirmed by the SF% values (see Tables S3 and S4, ESI†): ca. 28% (35%) for P–F and ca. 15% (25%) for P–H at the CC and MP2 (HF) levels. Moreover, the SF% difference between contributions coming from two directly bonded atoms was calculated passing from PX₅ to PX₃ allowing the variation of the polarisability of the bond to be evaluated; a decrease of ca. 16% was observed for P–F, whereas, for P–H a drastic decrease of ca. 3% —was observed (results based on CC and MP2 calculations; see Tables S8 and S9, ESI†).

HF results show an unbalanced polarisability from axial to equatorial bonds. QTAIM analysis relies on the knowledge of the charge density, which must be generated using some computational procedures, and if the charge density is built up using an insufficiently reliable method, its validity is limited. Thus, we do not expect HF-based results to be reliable, and only those results coming from CCSD and QCISD-based calculations will be discussed in the following.

The QTAIM atomic charges of PF₅ and PH₅ (see Tables S13 and S14, ESI†) show a high charge separation between P and X, with X = F (5.0 e) and H (3.5 e), which is a manifestation of high bond polarisation, as discussed by Cioslowski et al.⁴⁶ Ponec et al.²⁴ showed that such a kind of bond polarisation is typically associated with 3c-4e or 2c-2e highly polarised bonds, when the flexibility of the basis set is increased.

The evaluation of the dipole moment of F and H atomic basins, μ, in PF₅ and PH₅, (see Table S22, ESI†) shows for both H and F, axial positions to have smaller dipole moment than equatorial positions (viz., |μ| = 0.51 a.u. and 0.60 a.u. for F_ax and F_eq, respectively; |μ| = 0.26 a.u. and 0.36 a.u. for H_ax and H_eq, respectively). This confirms the higher bond polarisation of P–X_eq than P–X_ax.

To evaluate the average number of shared electron pairs, the δ(I,J) indexes were computed starting from HF calculations; in this case, HF theory gives the better results since the Fermi hole is the sole source of correlation between the electrons, as shown by Fradera et al. (1998).³⁷δ(I,J) indexes were compared with values that are normally expected for a single covalent bond, i.e., δ(I,J) ∼ 1.0,³⁷ and it was seen for both PH₅ and PF₅δ(I,J) < 1.0 (see Table S10, ESI†), confirming that there is a lower degree of electron sharing across these bonds. This is due to the highly-polarised nature of these bonds, especially for the two axial ones. A large bond polarisation, typical for a 3c-4e or two highly polarised 2c-2e bonds, depending on the basis set used,²⁴ clashes with the nature of a hydrogen atom, which preferably forms covalent bonds (characterised by low polarisation). Furthermore, the δ(I,J) values for the dissociation reactions for both PH₅ and PF₅ confirm an increase of electron sharing, which is larger for the former than the latter (Table S10, ESI†). This is due to the intrinsic electronic capability of hydrogen atoms to build shared shell interactions, which is in contrast with the partially ionic nature of the axial 3c-4e bond or two axial highly polarised 2c-2e bonds.²⁴ Ponec et al. (2004)²⁴ showed that increasing the basis set flexibility resulted in a significant change in the bonding scheme: 3c-4e transforms into a pattern of two more or less normal, albeit often very polar, two-center two-electron bonds.²⁴ Thus, it was decided to perform a DAFH analysis with small and large basis sets to better rationalise the results. A 3-centre index is a multicentre index; it therefore allows the presence of different types of bonds, such as 3c-4e or 3-centre-2-electron (3c-2e), to be identified depending on the magnitude of their calculated values in proximity to some theoretical limit proposed by Ponec et al. (1997).⁴⁷ On the other side, a 2-centre index allows the polarity of 2c-2e bonds to be determined going from 0.167 for H₂ (non-polar) to 0.483 for HCl (highest polar 2c-2e bond).^48,49 The enlarging of the basis set flexibility does not affect the bonding scheme of PH₅ which is always described by a 3 centre index as a 3c-4e bond for each basis set chosen. The negative value associated with H_ax–P–H_ax centres (3 centre index) equal to −0.011 (see Table 2) can explain its instability being far away from the limit for 3c-4e bonds (viz., −0.1875). On the other hand for PF₅ by enlarging the basis set flexibility the bonding scheme passes from 3c-4e to 2c-2e, highly polarised with a 2 centre index close to the limit proposed^48,49 (see Table 2). Furthermore, as recently shown by de Magalhães et al. (2016)⁵⁰ the values coming from 3 centre and 2 centre indexes can be used to understand which atoms will dissociate from a reagent during the reaction. In the PX₅ decomposition reaction the 3 centre and 2 centre indexes have values further away from the limit for X_ax than X_eq (see Table S23, ESI†) corroborating that the former, being more unstable, will dissociate from PX₅ to make X₂.

Table 2 Bader like calculated values of 3-centre bond indexes for X_ax–P–X_ax and 2-centre bond indexes for P–X_ax in PX₅

Bond indexes	Basis set
	6-31G**	6-311G**	cc-pVTZ
Hax–P–Hax	−0.013	−0.013	−0.011
Fax–P–Fax	0.010	—	—
P–Hax	0.654	0.643	0.636
P–Fax	0.341	0.434	0.323

---

F₂ has a quite low tabulated value of dissociation energy (159 kJ mol⁻¹)⁴⁹ due to the small atomic size and subsequently the small overlap of electron density of the two fluorine atoms. This increases the Pauli energy, E^Pauli, promoting a dipole–dipole interaction caused by the charge density rearrangement.⁵¹ This kind of interaction favours bonds that are not purely covalent. Due to its electronic nature, an F atom is therefore able to build a series of very highly polarised bonds with P in PF₅, wherein the electron sharing is reduced, since the expansion of the P valence shell plays a negligible role in bonding.⁵² Thus, the stability of the P–X_ax bond in PX₅ would imply a strong electronic correlation, in contrast with the weak electronic correlation of the hydrogen atom, which takes part in weakly shared shell interactions (viz., P–H_ax bonds). The shared shell interaction therefore causes instability in the PH₅ structure and thus dissociation is favoured by increasing electron sharing going to PH₃ and H₂ after dissociation.

Finally, to evaluate the stability by means of QTAIM, we tested the new topological descriptor introduced in this work – the Bader energy density, P_BADER (see eqn (2)). In the dissociation reactions of PX₅, P_BADER allowed both release/acquisition of energy associated with each atomic basin and also the volume variation to be taken into account. As explained by DAFH analysis for the PX₅ decomposition reaction (1); inside PX₅ the P–X_ax bonds are broken to give X₂, while P–X_eq bonds are distorted from the Cartesian xy-plane of the D_3h point group to the C_3V point group of PX₃. Thus, we calculated the P_BADER variation of an atomic basin like a change in its value passing from a reagent to a product: P_BADER(H_ax) and P_BADER(H_eq) show the same positive variation, ca. ΔP_BADER = 0.0006 a.u. per Bohr³; while P_BADER(F_ax) and P_BADER(F_eq) show opposite variation, ca. ΔP_BADER = −0.0133 a.u. per Bohr³ and ca. ΔP_BADER = 0.0610 a.u. per Bohr³, respectively (Table 3). The variation of atomic basin energy, ΔE, for F_ax and H_eq (ΔE = 0.63 a.u. and 0.29, respectively) does not suffice to explain the origin of PH₅ instability and PF₅ stability: in both cases we observe positive variation of atomic basin energy, which should be interpreted as a condition of instability.

Table 3 Volume (V), electronic energy (E) and Bader's energy density (P_BADER) of atomic basins

Atom	V (Bohr^3)	E (a.u.)	P BADER (a.u. per Bohr^3)	Atom	V (Bohr^3)	E (a.u.)	P BADER (a.u. per Bohr^3)
PX5
P	45.73	−339.47	−7.4239	P	14.15	−338.37	−23.9126
Heq	75.24	−0.86	−0.0114	Feq	106.95	−100.28	−0.9376
Hax	83.40	−0.88	−0.0106	Fax	106.95	−100.28	−0.9498
PX3
P	149.82	−340.23	−2.2710	P	108.67	−339.59	−3.1251
H	76.26	−0.82	−0.0108	F	114.31	−100.21	−0.8766
X2
H	59.02	−0.59	−0.0099	F	103.47	−99.65	−0.9630

---

A study of the Laplacian of charge density, as initially described by Bader²⁵ and later by Macchi,⁵³ can show the concentrated and depleted zones to describe the different types of bonding. For instance, due to high electron concentration close to the nucleus and the Fermi sea, transition metals may show an atomic basin through the plot of the null isosurface of the Laplacian of charge density. For other elements such a plot is not particularly informative and does not allow the visualisation of the atomic basins, showing only the border between concentrated and depleted zones. So a plot of the Laplacian of charge density for the PX₅ dissociation reactions by itself is unable to show the atomic basins, Table S24 (ESI†). In the present work, for the first time, thanks to the introduction of P_BADER, the volume of atomic basins assumes a primary character associated with the energy distribution inside of it. Focusing on volumes of different atomic basins of a P atom in different compounds studied here, it is possible to observe differences generated by ligands (Fig. 1). In PX₅ the atomic basin of the P atom is more “squeezed” when X = H than when X = F. In contrast, in PX₃ the atomic basin of the P atom does not show dramatic differences when X = H or F; the most remarkable difference consists of two extroflexions for H and only one for F.

Fig. 1 Plot of atomic basins of P with charge density depicted from highest concentration (white) to lowest concentration (black; 10⁻³ e Bohr⁻³).

Conclusions

A complete topological analysis of molecular systems involved in a PX₅ (X = H, F) dissociation reaction revealed profound differences in the bond type for P–F and P–H, classified as a transition zone and shared shell, respectively. The 3c-4e X_ax–P–X_ax bond avoids the P atom violating the octet rule, despite the shared shell interaction in the presence of a H atom, which is unfavourable for the electrostatic nature of this type of bond. The instability of PH₅ with respect to PF₅ was comprehensively shown thanks to the values of the new descriptor, P_BADER, which describes the release or acquisition of energy associated with each atomic basin along with the variation of the volume of the basin. P_BADER therefore represents a useful and intuitive quantity to describe the reaction course.

Conflicts of interest

There are no conflicts to declare.

Acknowledgements

The authors express their gratitude to Prof. Lukas Bucinsky and to Dr Gabriele Saleh for helpful discussions. This work was possible thanks to the computational facilities generously offered by Siberian Supercomputer Centre ICMMG SB RAS. Ch. T. would like to thank his supervisor Prof. Elena V. Boldyreva for ongoing support. E. B. expresses his gratitude to the Rector of Novosibirsk State University, and the CM of the RAS, Prof. M. P. Fedoruk, for his invitation to visit Novosibirsk State University.

References

1. C. Gatti and P. Macchi, Modern Charge-Density Analysis, Springer Dordrecht Heidelberg, London New York, 2012.
2. R. F. W. Bader, Atoms in molecules, A quantum theory, Oxford University Press, Oxford, 1990.
3. L. Pauling, The Nature of the Chemical Bond and the Structure of Molecules and Crystals: An Introduction to Modern Structural Chemistry, Cornell University Press, Ithaca, NY, 1960.
4. J. Yves, Molecular Orbital of Transition Metal Complexes, Oxford University Press, Paris, 2005.
5. H. Suzuki, Electronic Absorption Spectra and Geometry of Organic Molecules: An Application of Molecular Orbital Theory, Academic Press, NewYork, 1967.
6. R. D. Evans, The Atomic Nucleus, McGraw-Hill Book Company, New York, Toronto, London, 1955.
7. F. A. Carey and R. J. Sundberg, Advanced Organic Chemistry Part A: Structure and Mechanisms, Springer, Virginia, 2007.
8. W. A. Harrison, Electronic Structure and Properties of Solids: The Physics of the Chemical Bond, Dover Publication, New York, 1989.
9. C. Gatti and P. Macchi, Modern Charge Density Analysis, Springer, Dordrecht Heidelberg, London, New York, 2012.
10. G. N. Lewis, J. Am. Chem. Soc., 1916, 38, 762–785.
11. L. Pauling, J. Chem. Soc., 1948, 1461–1467.
12. R. Ponec and A. J. Duben, J. Comput. Chem., 1999, 20(8), 760–771.
13. C. Tantardini, E. V. Boldyreva and E. Benassi, J. Phys. Chem. A, 2016, 120, 10289–10296.
14. M. S. Schmøkel, S. Cenedese, J. Overgaard, M. R. V. Jørgensen, Y. Chen, C. Gatti, D. Stalke and B. B. Iversen, Inorg. Chem., 2012, 51, 8607–8616.
15. M. Durrant, Chem. Sci., 2015, 6, 6614–6623.
16. J. I. Musher, Angew. Chem., Int. Ed. Engl., 1969, 8, 54–68.
17. W. Kutzelnigg and J. Wasilewski, J. Am. Chem. Soc., 1982, 104, 953–960.
18. C. G. Pimintel, Physics, 1951, 19, 446; R. J. Hach and R. E. Rundle, J. Am. Chem. Soc., 1951, 73, 4321.
19. R. D. Hancourt, Lecture notes in chemistry, Springer, Berlin Heidelberg New York, 1982, vol. 30.
20. O. J. Curnow, J. Chem. Educ., 1998, 75, 10.
21. W. Kutzelnigg, Angew. Chem., Int. Ed. Engl., 1984, 23, 272.
22. C. A. Coulson, J. Chem. Soc., 1964, 1442.
23. Y. S. Cheung, Ch. Y. Ng, S. W. Chiu and W. K. Li, J. Mol. Struct.: THEOCHEM, 2003, 623, 1.
24. R. Ponec, G. Yuzhakov and D. L. Cooper, Theor. Chem. Acc., 2004, 112, 419–430.
25. R. F. W. Bader and C. Gatti, Chem. Phys. Lett., 1998, 287, 233.
26. C. Gatti and D. Lasi, Faraday Discuss., 2007, 135, 55.
27. E. Monza, C. Gatti, L. Lo Presti and E. Ortoleva, J. Phys. Chem. A, 2012, 115, 12864.
28. C. Gatti, Handbook Electron Density and Chemical Bonding, Springer, Berlin, Heidemberg, 2012.
29. M. S. Schmokel, S. Cenedese, J. Overgaard, M. R. V. Jorgensen, Y. S. Chen, C. Gatti, D. Stalke and B. B. Iversen, Inorg. Chem., 2012, 51, 8607.
30. L. Lo Presti, A. M. Orlando, L. Loconte, R. Destro, E. Ortoleva, R. Soave and C. Gatti, Cryst. Growth Des., 2014, 14, 4418.
31. F. J. Farrugia and P. Macchi, J. Phys. Chem. A, 2009, 135, 55.
32. C. Gatti, V. R. Saunders and C. Roetti, J. Chem. Phys., 1994, 101, 10686.
33. L. Bertini, F. Cragnoni and C. Gatti, Theor. Chem. Acc., 2007, 117, 847.
34. C. Gatti, F. Cragnoni and L. Bertini, J. Comput. Chem., 2003, 24, 422.
35. C. Gatti, P. Fantucci and G. Pacchioni, Theor. Chim. Acta, 1987, 72, 433.
36. E. Espinosa, I. Alkorta, J. Elguero and E. Molins, J. Chem. Phys., 2002, 117, 5529–5542.
37. X. Fradera, M. A. Austen and R. F. W. Bader, J. Phys. Chem. A, 1999, 103, 304–314.
38. R. F. W. Bader and Chérif F. Matta, J. Phys. Chem. A, 2004, 108, 8385–8394.
39. R. Ponec and D. L. Cooper, J. Phys. Chem. A, 2007, 111, 11294–11301.
40. F. W. BieglerKönig, R. F. W. Bader and T. Tang, J. Comb. Chem., 1982, 13, 317.
41. M. J. Frisch, G. W. Trucks, H. B. Schlegel, G. E. Scuseria, M. A. Robb, J. R. Cheeseman, G. Scalmani, V. Barone, B. Mennucci, G. A. Petersson, H. Nakatsuji, M. Caricato, X. Li, H. P. Hratchian, A. F. Izmaylov, J. Bloino, G. Zheng, J. L. Sonnenberg, M. Hada, M. Ehara, K. Toyota, R. Fukuda, J. Hasegawa, M. Ishida, T. Nakajima, Y. Honda, O. Kitao, H. Nakai, T. Vreven, J. A. Montgomery, Jr., J. E. Peralta, F. Ogliaro, M. Bearpark, J. J. Heyd, E. Brothers, K. N. Kudin, V. N. Staroverov, R. Kobayashi, J. Normand, K. Raghavachari, A. Rendell, J. C. Burant, S. S. Iyengar, J. Tomasi, M. Cossi, N. Rega, J. M. Millam, M. Klene, J. E. Knox, J. B. Cross, V. Bakken, C. Adamo, J. Jaramillo, R. Gomperts, R. E. Stratmann, O. Yazyev, A. J. Austin, R. Cammi, C. Pomelli, J. W. Ochterski, R. L. Martin, K. Morokuma, V. G. Zakrzewski, G. A. Voth, P. Salvador, J. J. Dannenberg, S. Dapprich, A. D. Daniels, Ö. Farkas, J. B. Foresman, J. V. Ortiz, J. Cioslowski and D. J. Fox, Gaussian 09, Revision D.01, Gaussian, Inc., Wallingford, CT, 2009.
42. AIMPAC. Available at: http://http:www.chemistry.mcmaster.ca/aimpac/imagemap/imagemap.htm.
43. T. Keith and R. F. W. Bader, Chem. Phys. Lett., 1993, 210, 223–231.
44. M. J. Frisch, G. W. Trucks, H. B. Schlegel, G. E. Scuseria, M. A. Robb, J. R. Cheeseman, J. A. Montgomery, Jr., T. Vreven, K. N. Kudin, J. C. Burant, J. M. Millam, S. S. Iyengar, J. Tomasi, V. Barone, B. Mennucci, M. Cossi, G. Scalmani, N. Rega, G. A. Petersson, H. Nakatsuji, M. Hada, M. Ehara, K. Toyota, R. Fukuda, J. Hasegawa, M. Ishida, T. Nakajima, Y. Honda, O. Kitao, H. Nakai, M. Klene, X. Li, J. E. Knox, H. P. Hratchian, J. B. Cross, V. Bakken, C. Adamo, J. Jaramillo, R. Gomperts, R. E. Stratmann, O. Yazyev, A. J. Austin, R. Cammi, C. Pomelli, J. W. Ochterski, P. Y. Ayala, K. Morokuma, G. A. Voth, P. Salvador, J. J. Dannenberg, V. G. Zakrzewski, S. Dapprich, A. D. Daniels, M. C. Strain, O. Farkas, D. K. Malick, A. D. Rabuck, K. Raghavachari, J. B. Foresman, J. V. Ortiz, Q. Cui, A. G. Baboul, S. Clifford, J. Cioslowski, B. B. Stefanov, G. Liu, A. Liashenko, P. Piskorz, I. Komaromi, R. L. Martin, D. J. Fox, T. Keith, M. A. Al-Laham, C. Y. Peng, A. Nanayakkara, M. Challacombe, P. M. W. Gill, B. Johnson, W. Chen, M. W. Wong, C. Gonzalez and J. A. Pople, Gaussian 09, Revision D.01, Gaussian, Inc., Wallingford CT, 2004.
45. F. Feixas, X. Girones, R. Ponec and J. Roithova, WinBader, v1.0, ICPF, AS CR, Prague, Czech Republic, 2003, available upon request, contact: PONEC@icpf.cas.cz.
46. J. Cioslowski and P. R. Surján, THEOCHEM, 1992, 255, 9–33.
47. R. Ponec and I. Mayer, J. Phys. Chem. A, 1997, 101, 1738–1741.
48. D. L. Cooper, R. Ponec, T. Thorsteinsson and G. Raos, Int. J. Quantum Chem., 1996, 57, 501.
49. J. J. De Corpo, R. P. Steiger, J. L. Franklin and J. L. Margrave, J. Chem. Phys., 1970, 53, 936.
50. Halua Pinto de Magalhães, Hans Peter Lüthi and Patrick Bultinick, Phys. Chem. Chem. Phys., 2016, 18, 846–856.
51. L. E. Forslund and N. Kaltsoyannis, New J. Chem., 2003, 27, 1108–1114.
52. M. Häser, J. Am. Chem. Soc., 1996, 118, 7311–7325.
53. P. Macchi and A. Sironi, Chemical bonding in transition metal carbonyl clusters: complementary analysis of theoretical and experimental electron densities, Coord. Chem. Rev., 2003, 238–239, 383–412; R. Ponec, Croat. Chem. Acta, 1994, 67, 55.

---

Footnote

† Electronic supplementary information (ESI) available: Computational details, molecular geometries, electronic energies, source function percentage values, delocalisation indexes, Espinosa indexes, and domain-averaged-Fermi hole analysis. See DOI: 10.1039/c7cp06130g

---