Electron Density Characteristics in Bond Critical Point (QTAIM) versus

---

ARTICLE pubs.acs.org/JPCA

Electron Density Characteristics in Bond Critical Point (QTAIM) versus Interaction Energy Components (SAPT): The Case of Charge-Assisted Hydrogen Bonding Barbara Bankiewicz,† Piotr Matczak,‡ and Marcin Palusiak*,‡ † ‡

Department of Theoretical Chemistry, University of Biazystok, Hurtowa 1, 15-399, Biazystok, Poland Department of Theoretical and Structural Chemistry, University of y
od
z, Tamka 12, 91-403 y
od
z, Poland

bS Supporting Information ABSTRACT: Charge-assisted hydrogen bonds (CAHBs) of N
H 3 3 3 Cl, N
H 3 3 3 Br, and P
H 3 3 3 Cl type were investigated using advanced computational approach (MP2/aug-ccpVTZ level of theory). The properties of electron density function defined in the framework of Quantum Theory of Atoms in Molecules (QTAIM) were estimated as a function of distance in H-bridges. Additionally, the interaction energy decomposition was performed for H-bonded complexes with different H-bond lengths using the Symmetry-Adapted Perturbation Theory (SAPT). In this way both QTAIM parameters and SAPT energy components could be expressed as a function of the same variable, that is, the distance in H-bridge. A detailed analysis of the changes in QTAIM and SAPT parameters due to the changes in H 3 3 3 A distance revealed that, over some ranges of H 3 3 3 A distances, electrostatic, inductive and dispersive components of the SAPT interaction energy show a linear correlation with the value of the electron density at H-BCP FBCP. The linear relation between the induction component, Eind, and FBCP confirms numerically the intuitive expectation that the FBCP reflects directly the effects connected with the sharing of electron density between interacting centers. These conclusions are important in view of charge density studies performed for crystals in which the distance between atoms results not only from effects connected with the interaction between atomic centers directly involved in bonding, but also from packing effects which may strongly influence the length of the H-bond.

’ INTRODUCTION Hydrogen bonding (H-bonding) plays an essential role in many physical, chemical and biological processes.1
8 Among different types of H-bonds those assisted by additional effects are the strongest. A good example can be the resonance-assisted hydrogen bonding (RAHB) that, according to the original concept proposed by Gilli et al.,9,10 is relatively stronger than other H-bonds because of the effect of resonance proceeding along the sequence of covalent bonds linking the protonaccepting group with the proton-donating one. In other words the contribution of charge-separated structures additionally stabilizes this type of H-bonding. Another example of the H-bonding strengthened by an additional effect is the so-called chargeassisted hydrogen bonding (CAHB).11
17 In this case the strengthening of H-bond results from the specific distribution of the charge in ion
ion or ion
neutral molecule complexes. Therefore, the different types of CAHBs can be defined on the basis of formal charge distribution in H-bridge (see Scheme 1). The CAHBs(+) are classified as interactions in which the cationic proton-donating fragment interacts with the formally neutral proton-accepting group. Consequently, the CAHBs(
) are the H-bridges in which the formally neutral proton-donating r 2011 American Chemical Society

group interacts with the anionic proton-accepting group. The most interesting situation appears when the cationic protondonating group interacts with the anionic proton-accepting group. The H-bonds of this type are classified as the so-called salt bridges, also known as CAHB(+/
) or double charge-assisted hydrogen bonds. Since CAHBs are stronger than nonchargeassisted H-bonds (non-CAHBs), they obviously have a larger impact on physical and chemical properties of the chemical species involved in complexation. For the strongest H-bonds the partially shared character of interaction was observed.18 This issue is addressed in detail in a recent review on covalency of H-bonding.19 Quantum Theory of Atoms in Molecules20,21 (QTAIM) gives the unique possibility to have an insight into a region of a system on the basis of physical properties of that system, i.e. it gives a possibility to divide the system (e.g., a molecule) into subsystems (e.g., atoms) on the basis of zero-flux in the electron density gradient field. Since the majority of chemists are interested in the Received: June 8, 2011 Revised: December 16, 2011 Published: December 17, 2011 452

dx.doi.org/10.1021/jp210940b | J. Phys. Chem. A 2012, 116, 452–459

The Journal of Physical Chemistry A

ARTICLE

studies on CAHBs78 it was shown that under the experimental conditions (e.g., in crystal structures) in CAHBs of different type and in their noncharge-assisted counterpart the H-bond length in H-bridge may vary over a large range of distances, from 1.8 to 4.0 Å. For this reason, the estimation of H-bonding strength directly from the experimental (X-ray) density in H-BCPs may in some cases be misleading. It was emphasized in many papers that the presence of a bond path (BP, the line linking the points of maximum electron density along the direction of the bond) linking a pair of atoms and the corresponding to it BCP fulfills the sufficient and necessary conditions that the atoms are bonded to one another (see, e.g., refs 80 and 81). This universal definition of any type of chemical bond is true when the forces acting on nuclei vanish, or in other words, when the system under consideration, which in specific conditions is represented by unique electron distribution, is in its equilibrium state. Since the properties of electron density function in BCPs are strongly related to the distance between atoms forming the bond, usually the direct comparison of such properties estimated for different bonds is not fully informative. For this reason, the representation of electron density in some specific point of space (e.g., in BCP) as a function of distance can be justified when one wants to have an insight into the changes in electron density function resulting from the change in distance between the interacting atoms. Such an analysis was already successfully applied to shared and closed-shell interactions.56
61,82 For a more detailed discussion on QTAIM-based parameters represented as a function of distance, see ref 83. In this paper we investigate different CAHBs of D—H 3 3 3 A (D = N, P and A = Cl, Br) type. Similarly to earlier paper,78 we estimate here the electron density properties at H-BCPs as a function of the distance. This gives us a unique possibility to compare in a conclusive way the QTAIM-based parameters of electron density function measured in BCP of CAHBs of different type. Additionally, we performed the decomposition of interaction energy, using the Symmetry-Adapted Perturbation Theory (SAPT)84,85 for the whole spectrum of distances in H-bridge, and tested the results against QTAIM data.

Scheme 1

relationship between the properties of some specific part of the molecule (e.g., reactivity center, substituent group, etc.) and the physical or chemical properties of that molecule (reactivity, accessibility to intermolecular interactions, etc.) the QTAIM theory becomes one of the most powerful tools in modern chemistry, forming a bridge between the advanced quantumchemistry and the experimental approach.22
25 The topological analysis of electron density can give valuable information about the properties of the system under consideration. For instance, it was found that electron density properties measured in the socalled bond critical point (BCP, the (3,
1) saddle point on electron density curvature being a minimum in the direction of the atomic interaction line and a maximum in the two directions perpendicular to it) may give information about the character of different chemical bonds, for example, covalent bonds,20,26,27 metal
ligand interactions,28
31 different kinds of H-bonds,4,32
36 including dihydrogen bonds37
41 and the so-called hydrogen
hydrogen interactions,42
47 as well as for other weak noncovalent interactions.48
51 It was also found that in the case of all those interactions there is a direct relation between the electron density at BCP and the energy of interaction between two centers forming the bond. For a small range of values this relation can be considered to be linear,4,52
55 although the true physical meaning has the exponential relation.56
59 Recently it was also shown that the changes in BCP because of the changes in distance are rather insensitive on the basis set and the method applied.60,61 It should be mentioned that the criticism of QTAIM was also formulated and several important cases in which QTAIM-based analysis may be misleading were also found.62
69 In general, the larger is the distance between two atoms forming the bond, the smaller is electron density in BCP of that bond and (usually) the weaker is the bond itself. For some specific cases QTAIM appeared to be useful in estimating the strength of H-bonding,4,19,33,58 particularly for intramolecular H-bonds for which the use of supermolecular approach is not possible.70
76 It was recently shown that the characteristics of electron density function in H-BCP are very strongly dependent on the type of centers (i.e., atoms or π-electrons) acting as proton acceptors in H 3 3 3 A bond (See for instance ref 52 and references therein and ref 77). On the other hand very recently the systematic studies78 on CAHBs of N
H 3 3 3 Cl type were reported and it was demonstrated there that the BCP characteristics change when the distance in H-bridge changes, but the location of formal charges in H-bridge does not practically affect these characteristics, despite the fact that the differences in interaction energies are dramatic, being larger than 120 kcal/mol as compared with CAHB(+/
) and its noncharge-assisted counterpart.78 In other words for the given (N)H 3 3 3 Cl distance the electron density and other parameters in BCP are practically the same, no matter whether non-CAHB or CAHB of any type, including CAHB(+/
), is considered. A similar conclusion was drawn for H
F bond on the basis of earlier studies on the charged and neutral HF molecular system.79 In the above-mentioned

’ METHODOLOGY All the structure calculations, were carried out with the Gaussian 09 suite of programs.86 Monomers, that is, D
H and A, and their complexes were optimized using the second-order Møller
Plesset perturbation theory (MP2)87,88 in conjunction with the Dunning’s correlation consistent basis set aug-ccpVTZ.89
91 For the CAHB(+) and CAHB(
) systems (see Scheme 1), no symmetry constraints were assumed in the course of full optimizations of these systems (i.e., the FOPT procedure in Gaussian 09). All the optimized CAHB(+) and CAHB(
) systems are in their stationary points that correspond to global minima on their potential surfaces. By contrast, the CAHB(+/
) complexes (see Scheme 1) were optimized with some restraints imposed on the D-H bond lengths, which prevented the initial geometries from converging into the H3D 3 3 3 H-A geometries. Thus, the obtained CAHB(+/
) systems, namely the D+— Hδ+ 3 3 3 A
type of H-bridges, neither exist as isolated molecular species, nor correspond to the lowest-energy structures. Instead, they characterize situations occurring commonly in experimental conditions, for example, ammonium and phosphonium halides crystals.92 Next we investigated the properties of all studied systems as functions of dH 3 3 3 A distance. Thus, we performed the scan along 453

dx.doi.org/10.1021/jp210940b |J. Phys. Chem. A 2012, 116, 452–459

The Journal of Physical Chemistry A

ARTICLE

the dH 3 3 3 A distance, keeping selected geometry parameters frozen. Namely, in the CAHB(+) and CAHB(
) systems, the D—H 3 3 3 A angle was kept frozen. Since it was necessary to prevent rotations of the interacting fragments which would lead to structures being far from the requested H-bridge, additional geometry constraints during the scan were enforced. QTAIM was used to analyze topological properties of the CAHBs in the investigated systems. The analysis made use of post-SCF MP2 densities (the keyword Density=Current in Gaussian 09). This analysis was performed with the help of the AIM200093,94 and AIMAll95 programs. The SAPT calculations were carried out using the aug-ccpVTZ basis set at the SAPT level equivalent to the supermolecular Many-Body Perturbation Theory (MBPT) through fourth order. These calculations were performed in SAPT2008.2 program.96 The SAPT formalism allowed us to decompose the interaction energy between the D
H and A monomers into a series of energy correction terms originating from various physical effects. In general, the SAPT interaction energy (ESAPT int ) is expanded in a perturbative series of the energy correction terms of increasing order. Treating the D—H 3 3 3 A complex in a perturbational manner assumes that the isolated D
H and A monomers are the zeroth-order approximation and the interaction between them is perturbing potential. In our case the SAPT expansion was truncated and took account of the energy correction terms up to the second order with respect to intermolecular interaction operator. ð10Þ

ð10Þ

ð20Þ

ð20Þ

¼ Eelst þ Eexch þ Eind, resp þ Eexch
ind, resp ESAPT int ð12Þ

ð22Þ

ð22Þ

ð20Þ

þ Eelst, resp þ t E ind þ t E exch
ind þ Edisp ð20Þ

ð13Þ

ð1Þ

þ Eexch
disp þ Eelst, resp þ εexch ðCCSDÞ ð2Þ

þ εdisp ð2Þ þ δEHF int, resp

is excluded from the discussion in the next sections. The assignment of the energy correction terms to four fundamental components of the interaction energy seems to be, to a certain extent, arbitrary and, therefore, a number of grouping schemes were proposed in the literature.97
100 In this work we adopted the grouping scheme (eqs 3
6) similar to that from ref 97. All further details concerning the SAPT part of calculations can be found in Supporting Information associated with this paper. The QTAIM parameters of H-BCPs and the SAPT decompositions were calculated within a wide range of distances between the H and A atoms involved in the D—H 3 3 3 A bonds. For each case of the H-bridge the range of distances taken into account was selected in the following way; the shortest dH 3 3 3 A distance was always an equilibrium distance found in optimization procedure, as explained above, whereas the maximum value of dH 3 3 3 A was the distance approximately close to that at which H-BCP vanishes (is no longer present). In this way, the minimum and maximum values of dH 3 3 3 A were estimated on the basis of properties being characteristic of each H-bridge under investigation. All dH 3 3 3 A distances are listed together with other parameters in Supporting Information Tables S1
S10.

’ RESULTS AND DISCUSSION To facilitate interpretation, we divide our discussion into two sections. In the first section we discuss the results for each type of CAHB separately, in the second section we provide more general observations and conclusions. In the both sections we will refer to the absolute values of the interaction energy and its components, keeping in mind that these energy parameters which correspond to stabilizing effects, have a negative value. All numerical data can be found in Supporting Information Tables S1
S10. In the main document we show only selected results using graphical representations. Influence of Atom Types in H-Bridge on the Energy and QTAIM Parameters of CAHBs in Question. The three types of

ð1Þ

CAHB(
) were investigated, namely, N—H 3 3 3 Cl
, P—H 3 3 3 Cl
, and N—H 3 3 3 Br
. In each case the halogen atom acts as an anionic proton-accepting center. From the comparison of interaction energies it can be said that Br
is a slightly weaker proton acceptor than Cl
. This is a generally known fact that was confirmed experimentally a long time ago.101
103 Similarly, the N
H fragment is a stronger proton donating center than the P
H bond. These observations are in agreement with general knowledge of H-bonding, which was recently formulated in the renewed definition of H-bond.104,105 If we compare the Eelst, Eind, and Edisp components of ESAPT int , it can be seen that for the N—H 3 3 3 Cl
and N—H 3 3 3 Br
systems the Eelst component is the largest in the whole range of dH 3 3 3 A. The Eelst component is followed by the Eind one, but the difference between them is always rather small (