The Concept of Density Functional Theory Based Descriptors and its Relation with the Reactivity of Molecular Systems: A Semi-Quantitative Study

Abstract

In this present paper, we have made an attempt to explain the theoretical basis for the empirical hardness/softness concepts to address the reactivity of molecular complexes in a semi-quantitative way within the framework of density functional theory. A model based on local hard-soft-acid-base principle has been proposed. The results obtained using some prototype charge transfer complexes, Lewis acid-base complexes and hydrogen-bonded complexes as examples, are in good agreement with the standard ab initio values. Although the model contains ad hoc parameters, it may form the basis of semi-quantitative description of inter-molecular interactions using hardness/softness parameters. The limitation, weakness and other critical issues of the present model are also discussed.

I. Introduction

One of the most important questions connected with the problem of reactivity of molecules in different environmental conditions is the prediction and interpretation of the preferred direction of a reaction and the product formation [1,2]. The study of molecular interactions has been a great challenge from the experimental and theoretical point of view [3]. There have been a lot of attempts to explain the nature of bonding and reactivity of molecular systems based on some intuitive ideas and empirical rules that are essentially derived from several experimental observations and many chemical facts [4]. During the development of the quantum chemical methods, many of the empirical chemical concepts were derived rigorously and it has provided a method for the calculation of the properties of chemical systems and the bonding that is involved in the molecular systems [4,5,6]. There also exist different kinds of theoretical approaches in correlating the reactivity of molecular systems based on different quantities, like molecular orbital density, charge on atoms, bond order, etc. [7]. In particular, the correlation between the frontier molecular orbitals and reactivity was brought out by the seminal works of Fukui and co-workers and it has been widely used to predict the nature of photochemical reactions. In the same period, Pearson had introduced the concept of hardness and softness (η/S) parameters in the context of explaining the reactivity of acids and bases [8]. Accordingly, he has systematized the reactivity of the acids and bases in terms of softness and hardness indices, based on the experimental observations. This study has eventually led to the celebrated principle, the so-called Pearson's Hard-Soft Acid-Base (HSAB) principle. It says that there is an extra stabilization when the soft acid combines with soft base and hard acid combines with hard base. HSAB principle has been very successful in rationalizing most of the acid-base reactions and it has become popular among the community of chemists because of its wider applicability and simplicity.

The nature of these basic chemical concepts, hardness and softness, called global reactivity descriptors (GRD), has been theoretically justified within the framework of density functional theory (DFT) [9,10]. Along with these global descriptors, other important local reactivity descriptors (LRD), such as Fukui function and local softness, were also proposed to rationalize the reactivity of a particular site in a molecule [11]. In the subsequent years, many groups have attempted to validate and prove the HSAB principle using the global and local reactivity descriptors [12]. These studies have led to some important insights about the nature of the reactivity and the stability of molecular systems in terms of η/S parameters [13,14,15,16,17,20,21,22,23,24,25]. Despite its partial success in describing the reactivity of the chemical systems, such studies have remained primarily qualitative. It should also be noted that the transformation of this qualitative principle into a quantitative form has been one of the challenging and difficult problems. This particular issue has been emphasized by many groups in the context of explaining the relative bond strengths of acid-base complexes or the reaction energies [6,18,19]. Drago and co-workers have given an expression to calculate the stabilization energy of the complexes in terms of the covalent and ionic (C and E) parameters [18a-b]. However, these are largely dependent on empirical parameters. In a recent study, we have raised some important questions that are pertinent for the establishment of a formal theoretical relationship between the structure and reactivity of molecular systems using the GRD and LRD [20,21]. In particular, we have attempted to provide a quantitative approach to calculate the interaction energy (IE) between the molecular systems using the density based descriptors. The quantitative approach is essentially based on the energy density perturbations up to the second order, as a functional of the number of electrons and the external potential and it has been formulated by Gazquez and Mendez [22,23] and Li and Evans [24]. Gazquez has calculated bond energies or IE for several diatomic molecular systems using GRD [16a]. He has also shown that the activation energy of a chemical reaction depends mainly on the difference between the hardness of the initial state of a reaction and the hardness of the transition state [16b]. Ghanty and Ghosh have also attempted to calculate the bond energies of diatomic molecular systems [15]. Their model involves the valence orbital radii, the electron density and the geometric mean of the homonuclear bond energies of the constituent atoms. The calculated bond energies for simple hetero-nuclear diatomic molecules are shown to agree very well with the experimental values. These models are essentially derived from the model proposed by Pauling that determines the bond energies in terms of the electronegativity [4]. All these models have been formulated to calculate the bond energy for some simple diatomic molecules in terms of the chemical potential and η/S parameters. For the case of poly-atomic molecular systems, the models are not directly applicable and it requires many parameters which are empirical in nature.

The model that is described in the present paper is to calculate the IE of the weak to moderate types of intermolecular interactions is based on local HSAB principle, developed by Gazquez and co-workers [22,23]. However, the formula contained an ad hoc parameter λ, which can not be rigorously calculated. In our recent work, we have made a critical study on the applicability of the local and global reactivity descriptors to describe weak interactions using the parameter λ as the charge transfer [20,21]. This has been tested for different weak interaction cases, for e.g. adsorption of molecules (N₂, CO₂ and CO) on cation exchanged zeolite surface [20]. It has also been applied to explain the energetics of nitrogen and oxygen molecules inside the zeolite cage and the results explain the selective sorption of N₂ molecule from air by the Na and Ca ion exchanged zeolite-A [25]. Although the model contains an ad hoc parameter and is thus semi-quantitative in nature, it has been successful in determining the interaction energies of these complexes. In view of this, a further study is yet to be made to understand the effect of the η/S parameters in the molecular interactions.

Accordingly, to study the above factors and to clarify the issues as detailed above, we have considered a wide range of systems that encompasses some of the chemical reactions, for instance, Lewis acid-base interaction, charge transfer interaction and hydrogen bonded (H-bonding) complexes. Another aspect of the present study is to demonstrate that the present model can be considered as a tool to monitor the influence of these descriptors in determining the reactivity of several types of complexes.

The paper is organized as follows: in section II the theoretical basis of the concepts of GRD and LRD, and the expression for the IE is given. Section III deals with the methodology and computational details and in Section IV, we discuss the results of the present study on the description of molecular interactions in terms of the η/S concepts.

II. Theoretical Background

II. 1 The global quantities

In DFT, the ground state energy of an atom or a molecule in terms of its electron density ρ(r) is written as [9,26],

E[ρ] = F[ρ] + ʃdrv(r)ρ(r)

(1)

where v(r) is the external potential that includes the nuclear potential also, and F[ρ] is the universal Hohenberg - Kohn functional composed of the electronic kinetic energy and the electron-electron repulsion energy. The first and second partial derivatives of E[ρ] with respect to the number of electron N under the constant external potential v(r) are defined as the chemical potential μ and the global hardness η of the system respectively [9,10].

$$\mu ={\left(\frac{\partial E}{\partial N}\right)}_{v(r)}$$

(2)

$$\eta =\frac{1}{2}{\left(\frac{{\partial }^{2}E}{\partial N^2}\right)}_{v(r)}$$

(3)

The inverse of the hardness is expressed as the global softness,

S= 1/2η

(4)

The global descriptor of hardness has been an indicator of overall stability of the system. It has been customary to use a finite difference approximation for μ and η. Using the energies of N, (N+1) and (N-1) electron systems, we get the operational definition of μ and η as [9],

μ ≈ - (IP + EA) / 2

(5)

η ≈ (IP – EA) / 2

(6)

where IP and EA are the first vertical ionization energy and electron affinity of the chemical species respectively.

II.2 Local Quantities:

The site-selectivity of a chemical system, can not, however, be studied using the global descriptors of reactivity. For this, appropriate local descriptors need be defined. An appropriate definition of local softness s(r) is given by [11],

$$s(r)={\left(\frac{\partial \rho (r)}{\partial \mu }\right)}_{v(r)}$$

(7)

such that,

$${\displaystyle \int s(r)dr=S}$$

Combining Eqs. 7 and the definition of global softness, we can write

$$s(\mathrm{r})={\left(\frac{\partial \rho (r)}{\partial \mu }\right)}_{v(r)}$$

(8a)

$$s(r)={\left(\frac{\partial \rho (r)}{\partial N}\right)}_{v(r)}{\left(\frac{\partial N}{\partial \mu }\right)}_{v(r)}$$

(8b)

= f(r) S

(8c)

where f(r) is defined as the Fukui function [11]. The Fukui function is defined as,

$$f(r)={\left(\frac{\partial \rho (r)}{\partial N}\right)}_v={\left(\frac{\partial \mu }{\partial v(r)}\right)}_N$$

(9)

The second relation can be obtained using the relation that density is the first partial derivative of energy with respect to external potential at constant N. The Fukui function describes the sensitivity of the chemical potential of a system to a local external potential. Using left and right derivatives with respect to the number of electrons, electrophilic and nucleophilic Fukui function and local softness can be defined. To describe site selectivity or reactivity of an atom in a molecule, it is necessary to condense the values of f(r) and s(r) around each atomic site into a single value that characterizes the atom in a molecule. This can be achieved by electronic population analysis. Thus for an atom k in a molecule, depending upon the type of electron transfer, we have three different types of condensed Fukui function of the atom k,

f_k+ = [q_k (N+1) - q_k (N)]  for nucleophilic attack

(10a)

f_k- = [q_k(N) - q_k (N-1)]   for electrophilic attack

(10b)

f_k0 = [q_k(N+1) - q_k(N-1)]  for radical attack

(10c)

where q_k is the gross electronic population of atom k in the molecule. The corresponding condensed local softnesses s_k+, s_k- and s_k0 can be defined. Parr and Yang proposed that larger value of Fukui function indicate more reactivity [11]. Hence greater the value of the condensed Fukui function, the more reactive is the particular atomic center in the molecule. Subsequently, Gazquez and Mendez proposed a local version of HSAB principle, which generally states that the interaction between any two chemical species will occur through the centers with nearly equal condensed Fukui function [22,23]. This can determine the behavior of different reactive sites with respect to the hard and soft reagents.

II.3 The Expression for the Interaction Energy

Using energy as a functional of number of electrons (N) and the external potential (v), the IE is defined as the difference of energies between the two interacting model systems A and B, and it can be defined in terms of two different steps, as [20,22,23],

ΔE_int = ΔE_v + ΔE_μ

(11)

$$\Delta E_v={\left(\frac{{({\mu }_A-{\mu }_B)}^2}{{\eta }_A+{\eta }_B}\right)}_v$$

(12)

$$\Delta E_{\mu }=-\frac{1}{2}{N}^2{({\eta }_{AB}-{\eta }_{AB}^{\ast })}_{\mu }$$

(13)

where η_AB and η*_AB the hardness of the complex at the equilibrium and at the isolated state, respectively. For the details of the mathematical part of derivation for the eq. 12-14, one can refer the work of Gazquez and Mendez [22,23] and by us [20]. Here the interaction between the system A and B is assumed to take place in two steps. In the first step, the interaction takes place at constant external potential through the equalization of chemical potential (ΔE_v). In the second step, A and B evolve toward the equilibrium state through changes in the electron density of the global system at constant chemical potential (ΔE_μ). The second step is actually a manifestation of principle of maximum hardness. One can relate the difference in the hardness terms present in the above eq.13 to the softness of system A and B with a proportionality constant (K) [20]. Thus, eq. 13 can be now rewritten in terms of the softness of the systems A and B as,

$$\Delta E_{\mu }=-\frac{1}{2}{N}_{AB}^{2}K\left(\frac{1}{(S_A+S_B)}\right)=-\lambda [1/(S_A+S_B)]/4$$

(14)

Herein, we introduce an ad-hoc term λ as the product of 2N² and the proportionality constant K. This parameter λ can be computed rigorously only through the softness of the molecular complexes. However, an ad-hoc value of λ as 0.5 and 1.0 has been used in the literature, which has been successful to a limited extent [22b,27]. In our earlier study, we have defined the parameter λ as the changes in the electron densities of the systems before and after the interaction process that will give the effective number of valence electrons that has participated in the interaction process [20]. Thus an expression for the term λ can be written as the difference of electron densities of the system A before and after the interaction,

$${\lambda }_A={\displaystyle \sum _{i=1}^M{\rho }_{Ai}^{eq}}-{\displaystyle \sum _{i=1}^M{\rho }_{Ai}^0}$$

(15)

Alternately, the term λ can be defined as the difference of electron densities for the system B,

$${\lambda }_B={\displaystyle \sum _{j=1}^N{\rho }_{Bj}^{eq}}-{\displaystyle \sum _{j=1}^N{\rho }_{Bj}^0}$$

(16)

where the first terms in eq. 15 and eq. 16 refer to the sum of the electron densities of the each atom in A and B in the molecule AB respectively and the second terms of these equations refer to electron densities of isolated systems A and B. For most practical cases, λ is almost same as the difference of electron density on the reacting atom of A or B. Recently, the above model has also been used by Chatterjee and co-workers to study the reactivity of several cationic sites in di-octahedral clays and the obtained results are in agreement with the experimental values [14a].

If the interaction between the systems occur through the kth atom of A with the lth atoms B, one can express the total IE from the local point of view, as [20],

$$\Delta E_{\mathrm{int}}^{AB}={\left({\mu }_A-{\mu }_B\right)}^2{\left(\frac{s_{Ak}{s}_{Bl}}{s_{Ak}+s_{Bl}}\right)}_v-\frac{\lambda }{4}{\left(\frac{1}{s_{Ak}+s_{Bl}}\right)}_{\mu }$$

(17)

where, μ_A and μ_B are the chemical potential of the A-B, respectively. The s_Ak and s_Bl refer to the condensed local softness of the atom k in the system A and l of the systems B, respectively.

III. Methodology and computational details

We have performed Möller-Plesset second order (MP2) [28] quantum chemical calculations to examine the reactivity of molecular systems. All the monomers and molecular complexes were optimized without any symmetry constraints, using MP2 level of theory through the standard split valence basis set, 6-31G(d,p). The restricted HF method has been used for the energy calculations of neutral and for the corresponding anionic and cationic systems, the restricted open shell HF method has been performed. The condensed Fukui function and local softness for each reactive atom were calculated via eq.10 using Mulliken population analysis [29]. The ab initio molecular orbital calculations were carried out using the GAMESS system of programs on an IRIX-6.2 silicon graphics workstation [30]. The parameter λ was calculated using eq.15 through the Mulliken population scheme. The reactive atoms that we have considered in our study for each type of interactions are, hydrogen, boron and halogen atoms in oxide-HX, BH₃-NH₃, NH₃-X₂ (electrophilic centers), respectively. Similarly, the nucleophilic reactive centers are oxygen atom in case of oxides, nitrogen atom in BH₃-NH₃ and NH₃-X₂ complexes. Accordingly, we have evaluated the values of the local softness s_k⁺ or s_k^- and f_k⁺ or f_k^- for electrophilic or nucleophilic centers in the molecules.

IV. Results and discussions

IV.1 Hydrogen bonding cases

The importance of H-bonds has been realized in different areas of chemistry. It plays a key role in determining the properties of chemical systems ranging from weak to strong bio molecular systems, crystal engineering and in molecular recognition [31]. These weak and strong H-bonds are classified based on the nature of proton acceptor and donor groups. In this present study, we have considered a novel H-bonding interaction between the ammonium oxide and phosphonium oxide with HX, where X=F, Cl and CN, and it involves the interaction of -N=O and -P=O groups with HX. These types of molecule exhibits stronger H-bonding than the usual cases and the IE lies in the energy range of -7.0 to -12 kcal/mol. In a recent study, Alkorta and Elguero have studied these systems using DFT and MP2 methods and they have correlated the geometry, energetics and electronic characteristics with electron density and Laplacian of the electron density [32]. We will now try to evaluate the IE of these complexes using the η/S parameters. The values of global and local reactivity descriptors are tabulated in Table 1. The value of chemical potential for HX is lower than that of oxides and it implies that the electron will flow from oxides to HX during the complex formation. The nucleophilicity (s_k^-) character of oxygen atom increases when the methyl group is substituted in place of hydrogen in NH₃O and it is found that s_k^- value of oxygen in PH₃O is less than the s_k^- of oxygen in NH₃O. In other words, (CH₃)₃NO is softer than NH₃O, and the analogous phosphonium oxide becomes hard. The electrophilicity (s_k⁺) value of hydrogen in HX decreases in the order of HF > HCl > HCN and the s_k⁺ value of hydrogen atom is considerably low in case of HCN. Based on these values, one can predict the reactivity of these molecular systems. For hydrogen donors, the reactivity order follows as HF > HCl > HCN and in case of hydrogen acceptors (oxides), the reactivity order is (CH₃)₃NO > NH₃O > PH₃O. Hence, the IE of oxide-HF complexes is expected to be more than the other complexes. Table 2 presents the computed IEs of all H-bonded complexes along with the available theoretical results and the stability order follow as (CH₃)₃NO > NH₃O > PH₃O. In case of NH₃O-HX complex, the predicted IE value is less than the reference value and in other cases, the values are in remarkably agreement with other theoretical values. Although the local softness of hydrogen in HCN is less than that of HCl,

Table 1. Global and local quantities of oxides, hydrogen halides, halogens and Lewis acids and bases. (Values in atomic units); Bolded atoms are the reactive atoms. Here the nucleophilic (Sk^-) atoms are O, N, P and S and the electrophilic (Sk⁺) atoms are H, F, Cl, Br, I and B.

System	Chemical potential	Global Softness	Fukui function	Local Softness
H 3NO	-0.088	2.158	0.645	2.158
H3PO	-0.140	1.799	0.414	0.745
Me3NO	-0.019	2.302	0.628	1.446
HF	-0.187	1.281	0.962	2.464
HCl	-0.153	1.613	0.686	2.313
HCN	-0.172	1.444	0.160	0.463
NH3	-0.095	1.699	0.478	0.813
F-F	-0.246	1.678	0.500	0.839
Cl-Cl	-0.211	2.253	0.500	1.127
Br-Br	-0.201	2.050	0.500	1.025
I-I	-0.158	2.147	0.500	1.073
H3B	-0.216	1.844	0.450	0.913
H3P	-0.109	1.845	0.588	1.084
H2O	-0.131	1.553	0.665	1.033
H2S	-0.108	1.860	0.768	1.430

Table 1. Global and local quantities of oxides, hydrogen halides, halogens and Lewis acids and bases. (Values in atomic units); Bolded atoms are the reactive atoms. Here the nucleophilic (Sk^-) atoms are O, N, P and S and the electrophilic (Sk⁺) atoms are H, F, Cl, Br, I and B.

System	Chemical potential	Global Softness	Fukui function	Local Softness
H 3NO	-0.088	2.158	0.645	2.158
H3PO	-0.140	1.799	0.414	0.745
Me3NO	-0.019	2.302	0.628	1.446
HF	-0.187	1.281	0.962	2.464
HCl	-0.153	1.613	0.686	2.313
HCN	-0.172	1.444	0.160	0.463
NH3	-0.095	1.699	0.478	0.813
F-F	-0.246	1.678	0.500	0.839
Cl-Cl	-0.211	2.253	0.500	1.127
Br-Br	-0.201	2.050	0.500	1.025
I-I	-0.158	2.147	0.500	1.073
H3B	-0.216	1.844	0.450	0.913
H3P	-0.109	1.845	0.588	1.084
H2O	-0.131	1.553	0.665	1.033
H2S	-0.108	1.860	0.768	1.430

Table 2. The interaction energies for the hydrogen bonded complexes (Energies in kcal/mol)^a

System	λ	ΔEv	ΔEμ	ΔEtot	ΔEbTE
H3NO----HF	0.0926	-4.01	-5.54	-9.55	-16.30
HCl	0.0746	-1.63	-4.68	-6.31	-10.30
HCN	0.0752	-0.89	-7.27	-8.15	-9.18
H3PO-----HF	0.050	-0.65	-3.97	-4.61	-9.79
HCl	0.0614	-0.05	-5.20	-5.25	-6.63
HCN	0.0488	-0.12	-7.83	-7.95	-6.85
(Me)3NO--HF	0.0521	-11.74	-3.05	-14.79	-16.92
HClc	0.0960	-7.02	-5.89	-12.91	-14.38
HCN	0.0598	-2.94	-5.59	-8.53	-8.42

^a The value of λ are given in atomic units. ^b ΔE_TE is the available theoretical IE values. See the ref. 32; ^c The geometry of (Me)₃NO-HCl has been taken from the ref. 32

Table 2. The interaction energies for the hydrogen bonded complexes (Energies in kcal/mol)^a

System	λ	ΔEv	ΔEμ	ΔEtot	ΔEbTE
H3NO----HF	0.0926	-4.01	-5.54	-9.55	-16.30
HCl	0.0746	-1.63	-4.68	-6.31	-10.30
HCN	0.0752	-0.89	-7.27	-8.15	-9.18
H3PO-----HF	0.050	-0.65	-3.97	-4.61	-9.79
HCl	0.0614	-0.05	-5.20	-5.25	-6.63
HCN	0.0488	-0.12	-7.83	-7.95	-6.85
(Me)3NO--HF	0.0521	-11.74	-3.05	-14.79	-16.92
HClc	0.0960	-7.02	-5.89	-12.91	-14.38
HCN	0.0598	-2.94	-5.59	-8.53	-8.42

^a The value of λ are given in atomic units. ^b ΔE_TE is the available theoretical IE values. See the ref. 32; ^c The geometry of (Me)₃NO-HCl has been taken from the ref. 32

it is interesting to note that the IE is close to the HCl complexes and it can be explained by the large value of the factor λ

IV.2 Charge transfer complexes and Lewis acid-base complexes

In this subsection, we will discuss the effect of the local and global descriptors on the most widely known cases like, charge transfer and Lewis Acid-base complexes. The study of electron donor-acceptor or charge transfer complexes is of great significance and Mulliken has theoretically formulated these types of interactions in terms of the covalent and ionic parts of wave function [33]. In the present study, we have considered the interaction between ammonia-halogen complexes (NH₃-X₂, where, X=F₂, Cl₂, Br₂ and I₂) and BH₃ with NH₃, PH₃, H₂O and H₂S. This type of complexes has been known for long time and several experimental and theoretical investigations have been performed to correlate the strength of the interaction with the charge transfer [34,35,36]. First, we present the results of our calculation in case of the NH₃-halogen complexes. The local and global values are given in the Table 1 and one can see that the Fukui function values are same for all halogen atoms. The local softness value for iodine is greater than that of other halogen atoms and hence, it expected that the IE value should be greater for the NH₃-I₂ system. We have presented the IE values for NH₃-X₂ systems in Table 3. While changing the X₂ system from fluorine to chlorine, the change in IE is 3kcal/mol and in case of NH₃-Br₂, IE changes only one kcal/mol compared to NH₃-Cl₂ complex. On the other hand, for the NH₃-I₂ complex, IE is observed to be high and it could be due to the large value of λ. In general, it can be seen that the stability of the ammonia complex formed with X₂ increases with the increase of the atomic number or increasing its electronegativity difference. These results are consistent with the recent literature works and with other reactivity correlations [31,32,33,34,35,36].

We will now discuss the applicability of our methods to the case of Lewis Acid-Base Complexes (LABC) and herein, we have studied selective important complexes only for which other theoretical results are available. The values of the reactivity descriptors and the IE are tabulated in Table 1 and Table 4, respectively. It can be seen from the Table 1, that the local softness and Fukui function values is increased substantially with the replacement of nitrogen by phosphorous in ammonia and with oxygen by sulfur in water. Although the local softness value of phosphorous in PH₃ is high, the interaction of BH₃ with NH₃ is more favored by 5 kcal/mol due to the high value of the parameter λ for

Table 3. The interaction energies for the Charge transfer complexes (Energies in kcal/mol)^a

System	λa	ΔEv	ΔEμ	ΔEtot
NH3-F2	0.0163	-5.91	-1.55	-7.46
Cl2	0.0889	-4.01	-7.19	-11.19
Br2	0.1585	-3.21	-9.26	-12.46
I2	0.2079	-1.15	-17.29	-18.44

^a The value of λ are given in atomic units.

Table 3. The interaction energies for the Charge transfer complexes (Energies in kcal/mol)^a

System	λa	ΔEv	ΔEμ	ΔEtot
NH3-F2	0.0163	-5.91	-1.55	-7.46
Cl2	0.0889	-4.01	-7.19	-11.19
Br2	0.1585	-3.21	-9.26	-12.46
I2	0.2079	-1.15	-17.29	-18.44

^a The value of λ are given in atomic units.

Table 4. The interaction energies for the Lewis acid-base complexes (Energies in kcal/mol)

System	λa	ΔEμ	ΔEv	ΔEtot	ΔEbTE
BH3-NH3	0.3067	-3.94	-27.88	-31.82	-31.30
PH3	0.2729	-3.51	-21.43	-24.94	-20.56
H2O	0.2203	-2.19	-17.76	-19.95	-10.30
H2S	0.3021	-4.05	-20.23	-24.28	-12.02

^a The value of λ are given in atomic units.

^b ΔE_TE is the available theoretical IE values. See the ref. 36

Table 4. The interaction energies for the Lewis acid-base complexes (Energies in kcal/mol)

System	λa	ΔEμ	ΔEv	ΔEtot	ΔEbTE
BH3-NH3	0.3067	-3.94	-27.88	-31.82	-31.30
PH3	0.2729	-3.51	-21.43	-24.94	-20.56
H2O	0.2203	-2.19	-17.76	-19.95	-10.30
H2S	0.3021	-4.05	-20.23	-24.28	-12.02

^a The value of λ are given in atomic units.

^b ΔE_TE is the available theoretical IE values. See the ref. 36

BH₃-NH₃ than BH₃-PH₃ complex. In case of H₂S, the value of local softness of sulfur and the λ parameter are higher than the local softness of oxygen in H₂O and hence, the interaction of BH₃ with H₂S is more favored than H₂O.

In order to test the efficiency of our method, the IE calculated by our method is compared with the IE values obtained by the standard density functionals and MP2 method [35] and all the values are presented in Table 5a and Table 5b. It can be seen that VWN functional overestimates the IE (~30kcal/mol) in comparison to the other methods and the inclusion Hartree-Fock exchange term reduces the IE term dramatically. It is pertinent to note that the predicted IE values from our method are in agreement with the MP2 results. As the 6-31G(d,p) basis set is not available for iodine cases, we have employed 3-21G(d,p) basis set and it could be the reason that our model predicts high IE value. Similarly, in case of BH₃-NH₃ adduct case, the IE is compared with the experimental and other available highly accurate theoretical results [36]. There is an excellent agreement with the value obtained through our model and G2(MP2) methods and they differ 0.5 kcal/mol from the experimental result.

Table 5a. Comparison of the IE (in kcal/mol) values of NH₃-X₂ complexes obtained by our method and by other standard DFT and MP2 methods.^a

System	VWN	BP	PP	PW	HHLYP	MP2	Our method
NH3-F2	-35.8	-26.9	-33.4	-26.0	-2.4	-2.2	-7.47
Cl2	-20.1	-13.0	14.1	-13.1	-9.2	-7.8	-11.19
Br2	-21.3	-13.2	-14.5	-14.7	-10.1	-11.0	-12.46
I2	-18.0	-13.0	-13.7	-11.3	-10.2	-11.6	-18.44

^a All the abbreviations for DFT functionals can be referred in the ref. 35

Table 5a. Comparison of the IE (in kcal/mol) values of NH₃-X₂ complexes obtained by our method and by other standard DFT and MP2 methods.^a

System	VWN	BP	PP	PW	HHLYP	MP2	Our method
NH3-F2	-35.8	-26.9	-33.4	-26.0	-2.4	-2.2	-7.47
Cl2	-20.1	-13.0	14.1	-13.1	-9.2	-7.8	-11.19
Br2	-21.3	-13.2	-14.5	-14.7	-10.1	-11.0	-12.46
I2	-18.0	-13.0	-13.7	-11.3	-10.2	-11.6	-18.44

^a All the abbreviations for DFT functionals can be referred in the ref. 35

Table 5b. Comparison of the IE (in kcal/mol) of BH₃-NH₃ complex obtained by our method and by other high-level accurate quantum chemical calculations.^a

Theoretical methods	IE
G2(MP2)	-31.3
G2	-31.4
CCSDT	-31.1
NL SCF	-32.1
DFT	-32.8
MP2/TZP	-33.7
Our method	-31.8
Experiment	-31.1

^a All the abbreviations for the other theoretical methods can be referred in the ref. 36

Table 5b. Comparison of the IE (in kcal/mol) of BH₃-NH₃ complex obtained by our method and by other high-level accurate quantum chemical calculations.^a

Theoretical methods	IE
G2(MP2)	-31.3
G2	-31.4
CCSDT	-31.1
NL SCF	-32.1
DFT	-32.8
MP2/TZP	-33.7
Our method	-31.8
Experiment	-31.1

^a All the abbreviations for the other theoretical methods can be referred in the ref. 36

IV.3 The importance of chemical potential equalization process and the principle of maximum hardness

As we have detailed in section II, we will now examine the essential role of the chemical potential equalization and evolution of hardness process in the molecular interaction process. The analysis of the contribution of each process (ΔEv and ΔE_μ) in stabilizing the molecular complexes, reveals that the complexes are mostly driven by the process of the reshuffling or dynamical relaxation of the charge distribution at constant chemical potential. The only exceptional case is the (CH₃)₃NO-HF complex for which the complex is stabilized through the first term by an amount of -11.7 kcal/mol. One should also note that the second term ΔEv is largely dependent on the factor λ, the effective number of electrons that has been transferred from the reacting systems, and it is interesting to find some relationship between the factor λ and the soft-soft and hard-hard interactions. In general, the value of λ in soft-soft type of interactions is significantly higher than that hard-hard interaction. For instance, in case of NH₃-I₂, BH₃-NH₃, NH₃NO-HF (soft-soft), the value of λ is higher than the corresponding, NH₃-F₂, NH₃-PH₃, NH₃NO-HCl (soft-hard and hard-hard) complexes. The present analysis in relating the factor λ with different types of interactions, is quite analogous to Klopman's chemical reactivity theory [6]. Klopman has shown that the soft-soft interactions are mostly controlled by the frontier orbitals of the interacting systems and in case of the hard-hard type of interactions, the contribution of electrostatic interaction becomes more significant. Although, the present method and Klopman's chemical reactivity theory has different origin, the arguments are quite similar.

We have also carried out a study of basis set effects on the global and local reactivity descriptors by introducing the polarization and diffusion functions. We have checked that the IE evaluated using the present model is consistent with the value of IE calculated by the conventional quantum chemical methods in different basis sets. The details of this study can be referred elsewhere [37].

Before we conclude this section, we would like to mention the limitations of our present approach. In particular, the effectiveness and accuracy of the present method lie on the computation of the local descriptors and these are highly dependent on the basis set and level of theory that is used in the calculation. However, it should be noted that these issues are quite common in any kind of models and the accuracy will be ultimately determined by the level of computation. Despite the arbitrary nature of the population analysis and the basis set, that has been applied in the computation of each term present in our approach, we could still get the reliable IE values that are in agreement with the experimental or available theoretical results. These values can be further improved by making judicious choice of basis set and population methods. Recently, Roy et al has shown that Hirshfeld population scheme would be better choice for the computation of condensed Fukui function values and this method seems to be less sensitive to the basis set [38]. The problem of defining the factor λ is still an issue. Although the present definition of λ works well in most cases, it may not work for cases where the charge transfer becomes negligible or zero.

V. Conclusions

"At the present time, we do not know exactly how to calculate the total changes in energy due to the moving about of nuclei in a chemical reaction, using changes in chemical potential and hardness parameters as the determinates. But this does not appear to be an intractable problem. In the future, we can hope to have a much better understanding of the relationship between total energy changes and changes in chemical hardness. This should lead to an improved, more quantitative HSAB principle" [39]. This provocative comment made by Pearson illuminated an issue that was believed to represent a challenging problem. However, the present and earlier studies show that it can provide a reliable semi-quantitative model to calculate the IE between the molecular systems in terms of the η/S, chemical potential parameters along with the Fukui function (GRD and LRD) [20,21,25]. Although the present model has been successful in describing several weak to moderate types of molecular interactions, it involves an ad hoc parameter λ, which can not be computed rigorously. A further study should be made on the parameter λ and other parameters involved in the expression of IE. Future work should focus on the general classification of the type of interaction that is involved in a large number of complexes based on the mean values of these descriptors.