Comparison of the DFT-SAPT and Canonical EDA Schemes for the

Article pubs.acs.org/JCTC

Cite This: J. Chem. Theory Comput. XXXX, XXX, XXX−XXX

Comparison of the DFT-SAPT and Canonical EDA Schemes for the Energy Decomposition of Various Types of Noncovalent Interactions Olga A. Stasyuk,† Robert Sedlak,†,‡ Ceĺ ia Fonseca Guerra,§,∥ and Pavel Hobza*,†,‡ †

Institute of Organic Chemistry and Biochemistry, Academy of Sciences of the Czech Republic, 166 10 Prague 6, Czech Republic Regional Centre of Advanced Technologies and Materials, Department of Physical Chemistry, Palacký University, 771 46 Olomouc, Czech Republic § Department of Theoretical Chemistry and Amsterdam Center for Multiscale Modeling, VU Amsterdam, De Boelelaan 1083, 1081 HV Amsterdam, The Netherlands ∥ Leiden Institute of Chemistry, Gorlaeus Laboratories, Leiden University, P.O. Box 9502, 2300 RA Leiden, The Netherlands

Downloaded via NAGOYA UNIV on June 21, 2018 at 16:50:55 (UTC). See https://pubs.acs.org/sharingguidelines for options on how to legitimately share published articles.

‡

S Supporting Information *

ABSTRACT: Interaction energies computed with density functional theory can be divided into physically meaningful components by symmetry-adapted perturbation theory (DFT-SAPT) or the canonical energy decomposition analysis (EDA). In this work, the decomposition results obtained by these schemes were compared for more than 200 hydrogen-, halogen-, and pnicogen-bonded, dispersion-bound, and mixed complexes to investigate their similarity in the evaluation of the nature of noncovalent interactions. BLYP functional with D3(BJ) correction was used for the EDA scheme, whereas asymptotically corrected PBE0 functional for DFTSAPT provided some of the best combinations for description of noncovalent interactions. Both schemes provide similar results concerning total interaction energies and insight into the individual energy components. For most complexes, the dominant energetic term was identified equally by both decomposition schemes. Because the canonical EDA is computationally less demanding than the DFT-SAPT, the former can be especially used in cases where the systems investigated are very large.

■

simplified force fields or specific potentials. Studies of model complexes, such as water dimer, benzene dimer, etc. showed how reliable the utilized decompositions are, or which of the attractive forces is responsible for aggregation.4,9−13 The different variants of the decomposition schemes usually apply slightly different definitions for the energy components of the interaction energy. Moreover, separation of the charge transfer energy from the polarization part for short intermolecular distance is problematic, especially in the case when the extended basis set is utilized.13 The possible applicability of the energy decomposition is mostly limited by the fact at which level of theory particular scheme is implemented. The schemes based on the SAPT, with their coupled cluster implementation, can be utilized only for molecular systems consisting of several atoms up to couple dozens of atoms. In turn, combination of SAPT with DFT (DFT-SAPT), in which the DFT monomer description accounts for intramolecular correlation effects, can be applied for larger systems.14,15 To extend the applicability of DFTSAPT to even larger complexes, empirical dispersion corrections instead of computationally expensive dispersion

INTRODUCTION Importance of the noncovalent interactions in natural sciences is indisputable. They are responsible for many important processes in various fields, such as molecular chemistry, medicinal chemistry, biochemistry, and catalysis. Therefore, the interest in their investigation is still ongoing. Experimental techniques do not provide information about the driving force of the interactions. Moreover, there is no quantum chemical operator which would directly lead to the value of the interaction energy or at least to any of its components. However, it has been shown that energy decomposition analysis (EDA) can provide us this insight from computations. Different energy decomposition schemes have been developed to understand the chemical bonding.1−6 These methods build upon the original works by Morokuma7 (variational approach) or Jeziorski et al.8 (perturbative approach called Symmetry Adapted Perturbation Theory (SAPT)). It should be noted that in SAPT-based schemes the total interaction energy is constructed from different energy terms, whereas in EDA schemes the interaction energy is decomposed into energy contributions. The strength of all decomposition schemes is that they provide insight into the nature of the chemical interactions. Additionally, energy decomposition can also be helpful when designing new or © XXXX American Chemical Society

Received: January 12, 2018

A

DOI: 10.1021/acs.jctc.8b00034 J. Chem. Theory Comput. XXXX, XXX, XXX−XXX

Journal of Chemical Theory and Computation

■

and exchange-dispersion terms have been proposed.16,17 The field of EDA schemes is also rapidly growing, and new approaches which move applications toward much larger complexes are continuously developed. The investigations of small biomolecular clusters18,19 as well as larger complexes20 which are treated via QM/MM approach have already appeared in the literature. Some EDA schemes work within the framework of the fragment molecular orbital method.3,9,21,22 This approach divides the complex into several fragments which are treated independently. Because of this, it is possible to apply this method for large molecular clusters, such as proteins. Moreover, the periodic implementation of the EDA has been already introduced.23 In several papers a comparison of different decomposition methods has been provided.24−28 Description and comparison of various variational and perturbation based decomposition schemes have been reported by Phipps et al. in the comprehensive review for prototypical protein−drug interaction patterns.27 They noted that KM EDA,29 RVS EDA,30 and SAPT provides overall reasonably fair energy components, whereas the NEDA scheme31−33 yields sometimes extremely large and chemically unrealistic results for their systems. In turn, the ALMO EDA scheme9 was recognized as the best one for the systems relevant to drug optimization. Test of effective fragment potential (EFP) energy components and SAPT energies for the S22 data set showed that the main source of errors in EFP comes from Coulomb and polarization terms and requires further improvements.26 Introduced very recently, Local Energy Decomposition (LED)28 based on DLPNO− CCSD(T) method showed good agreement with DFT-SAPT for the total interaction energy, whereas the methods differ in their estimation of the dispersion contributions with systematically larger values for DFT-SAPT than for LED. In this study, we apply methods utilizing symmetry adapted perturbation theory (SAPT) in combination with density functional theory (DFT), namely, DFT-SAPT, and the canonical EDA as implemented in ADF software, from now on referred to as EDA. The first method is representative of the hybrid perturbative and supermolecular approach, whereas the second one of the variational approach. Both methods divide the total interaction energy into physically meaningful components which allow understanding the origin of interaction but using different approaches. In order to compare the effectiveness of particular decomposition schemes, we test both methods on a large set of noncovalently bound complexes (more than 200) and demonstrate for them the consistency of energy components and stress their distinguishing features. As model complexes, we chose dimers with hydrogen, halogen, pnicogen bonds, π···π stacking and CH···π contact, as well as with mixed type of interaction to cover various interaction types according to their nature (electrostatic, charge-transfer, dispersion or their combination) and wide range of stabilization energy. The aim of this study is not to provide rigorous comparison of the physical background of both studied methods but rather to compare their applicability to different types of noncovalent interactions (see the General Remarks section). The crucial feature of this study is the robustness of the comparison performed for the large set of data which includes various complexes with noncovalent interactions.

Article

METHODOLOGY

Both presented methods are implemented in two different program packages: DFT-SAPT in the Molpro34 and the canonical EDA in the ADF program. 35 The detailed description of both methods can be found in the original papers on DFT-SAPT36−38 and EDA22,39,40 and more recent reviews.6,41,42 DFT-SAPT. The original many-body SAPT method8 is computationally expensive. An application of the DFT method for accounting for the electron correlations within each monomer allows one to reduce significantly the computational cost. The DFT-SAPT method36−38 belongs to so-called “hybrid” methods as it combines both perturbation based (the framework in which the first- and second-order interaction energy terms are defined) and supermolecular approach (δHF term). There is also a similar approach, known as SAPT(DFT).43 The DFT-SAPT method decomposes quantitatively the interaction energy into several components of the first and second perturbation order and supermolecular δ(HF) term: (1) (1) (2) (2) (2) (2) + Eexch + Eind + Eexch Eint = Eelst − ind + Edisp + Eexch − disp

+ δ(HF )

(1)

Some of these terms can be grouped in order to indicate commonly understood physical quantities. The first term, E(1) elst , corresponds to the electrostatic interaction between monomers, E(1) exch contains the exchange−repulsion contribution, the (2) next two terms E(2) ind and Eexch−ind can be combined into the one named as an induction part, which covers also a charge-transfer (2) energy, and the sum (E(2) disp + Eexch−disp) accounts for dispersion contribution. The last term, δ(HF), indicates higher-order induction and exchange−induction terms in the Hartree−Fock energy and is determined as the difference between the HF interaction energy and the sum of all the first- and secondorder contributions (obtained with the HF wave functions), except for the dispersion and exchange−dispersion terms. This term can be grouped with the induction part.44 The importance of the δ(HF) correction for the accuracy of the total interaction energy is indisputable, especially for complexes with short intermolecular distances. The DFT treatment of monomers is performed with the asymptotically corrected exchange-correlation functional. The choice of functional is crucial as it also affects the results of couple perturbed response equations. It was previously shown that among DFT functionals the PBE0AC functional provides one of the best results for noncovalent interactions in comparison with coupled cluster based SAPT.45 In this work, the DFT part was treated using the PBE0 or LPBE0 functional with aug-cc-pVDZ or aug-cc-pVTZ basis sets as implemented in Molpro 2010 package.34 Heavy atoms, such as Br, Cl, As, Sb, and Bi, were treated by the pseudopotentials to cover the relativistic effects of the inner-core electrons. As it is known that DFT canonical virtual orbital energies are less than optimal for use in a perturbative scheme, their values were corrected before the SAPT treatment by a gradient-controlled shift procedure.46 This step uses the difference between the exact vertical ionization potential (IP) and the energy of the highest occupied molecular orbital (HOMO) of each monomer obtained by additional calculations of a neutral system and its cation. The DFT (gradient shift) calculations were done utilizing the Turbomole package.47 B


Article

Journal of Chemical Theory and Computation Canonical EDA. The intermolecular interaction energy ΔEint between two fragments is calculated by the supermolecular approach; i.e., the energies of the monomers EA and EB are subtracted from the energy of the dimer, EAB, as ΔEint = ΔEAB − ΔEA − ΔEB

perturbation expansion. This is quite obvious, assuming that the perturbation is within the SAPT theory defined as the intermolecular Hamiltonian, while the unperturbed Hamiltonian is defined as a sum of the Hamiltonians of isolated monomers. Demonstration of comparability of these two energy decomposition schemes is important because DFT-SAPT approach is computationally more expensive and can be applied to medium-sized systems, whereas canonical EDA is faster and allows calculating energy terms for large biological complexes with few hundreds of atoms. Moreover, EDA can be used for catalysis, reactions, formation of electron pairs, transition metals, etc. Finally, let us make a few comments about possible sources of discrepancies between the two energy decomposition methods. In this work we used two functionals (PBE0AC in DFT-SAPT vs BLYP-D3(BJ) in EDA), which were chosen as most commonly used functionals in these methods for noncovalently bound systems. Keeping in mind that the main purpose of the paper is qualitative comparison of the decomposition methods, we checked the correlations between EDA results for PBE0-D3 and BLYP-D3 functionals and DFTSAPT results for 10 representative complexes. Additionally, we checked the effect of geometry optimization on energetic terms. Description of complexes and obtained correlations are given in the Supporting Information. The BLYP-D3 functional provides very similar correlations for the particular energy components as PBE0-D3 and can be used without additional geometry optimization. Comparison of the Energy Components: DFT-SAPT vs Canonical EDA. The direct comparison of specific energy terms is difficult as they are defined in different way. Since the (1) electrostatic, E(1) elst and ΔVelstat, and the repulsion, Eexch and ΔEPauli, components are defined within both methods in similar way, their direct comparison can be done without any modifications. Unfortunately, this is not the case for the rest of the components. (2) Considering the dispersion components, Edisp=(Edisp + ) from DFT-SAPT and ΔE from the canonical E(2) exch−disp disp EDA, it is evident that their direct comparison is not justified. On one hand, within the DFT-SAPT method, the dispersion term is represented as a solution of the frequency-dependent coupled-perturbed Kohn−Sham equations and describes all intermolecular correlation effects to the second order. On the other hand, within the canonical EDA method, it is mostly long-range atom−atom pair empirical correction introduced by Grimme. Regarding magnitudes of the dispersion component, one can expect larger values for the DFT-SAPT dispersion contrary to the EDA one, particularly in the case when a significant overlap of the monomer electron densities takes place. This is due to the fact that the Grimme’s correction almost exclusively covers only the long-range correction to the dispersion and does not include its short- and middle-range parts, which are contained in GGA functional. Furthermore, D3 dispersion correction was developed at the def2-QZVP basis set to provide results quite close to the KS limit; therefore, this correction can be directly applicable without any scaling. However, DFT-SAPT dispersion energy converges more slowly with basis set size than the remaining parts of the interaction energy but can be successfully extrapolated to the CBS limit applying a scaling factor f.61 This factor can be obtained by comparing the dispersion energy obtained at each basis set with CBS limit and can be

(2) 22,39,40

According to the Ziegler−Rauk method, the interaction energy, ΔEint, can be decomposed into three physically meaningful components representing different steps toward the formation of a complex from two selected fragments: ΔEint = ΔVelstat + ΔEPauli + ΔEoi

(3)

The first term, ΔVelstat, corresponds to the classical electrostatic interaction between the fragments superimposed with their frozen charge distribution (ρA and ρB) at the geometry of the complex and is usually attractive. The Pauli repulsion energy, ΔEPauli, comprises the destabilizing interactions between occupied orbitals and is responsible for the steric repulsion originating from the Pauli antisymmetry principle. The orbital interaction energy, ΔEoi, accounts for the charge transfer (donor−acceptor interactions between the occupied molecular orbitals on one fragment with the unoccupied molecular orbitals of the other fragment) and intrafragment polarization (mixing of occupied and virtual orbitals within the same fragment). These polarization and donor−acceptor interactions are difficult to separate as fragment orbitals can be involved in both interactions. The canonical EDA has been further developed by Bickelhaupt and Baerends for the analysis of electron-pair bonding (formation of doubly occupied molecular orbital (MO) from singly occupied molecular orbitals) and other types of bonding involving open-shell fragments,48,49 but this is beyond the scope of this study. Later, the empirical term ΔEdisp was added to the interaction energy to account for the long-range dispersion effects, as introduced by Grimme and co-workers.50,51 This term is computed independently from the Kohn−Sham energy to yield the total interaction energy of the system. All calculations were carried out using the Amsterdam Density Functional (ADF) program.35 The interaction energies were calculated using the generalized gradient approximation (GGA) with BLYP functional52,53 and DFTD3(BJ) correction51 as the best correction for noncovalent interactions.54,55 Combination BLYP-D3(BJ) was recently shown as one of the best DFT methods for noncovalent interactions.56,57 The MOs were expanded in uncontracted sets of Slater type orbitals (STOs) containing diffuse functions with two sets of polarization functions (TZ2P).58 For halogen and pnicogen bonded complexes, relativistic correction were introduced by zero-order regular approximation (ZORA).59,60 General Remarks. Let us make a short comment about the limitations of particular schemes with respect to the intermolecular distances. As already mentioned, both methods are commonly tested for the noncovalent interactions for which they provide relatively reasonable results in the vast majority of cases. The canonical EDA scheme provides most likely reasonable results regardless the intermolecular distance. On the other hand, it is known that the perturbation based SAPT method sometimes fails to describe accurately noncovalent interactions with short intermolecular distances.16 Explanation of this drawback can be non-negligible size of perturbation with respect to the unperturbed Hamiltonian, which is a necessary condition for the convergence of the C


Article


Table 1. Errorsa of the Methods with Respect to the Benchmark CCSD(T)/CBSb Results for Different Types of Interactionsc canonical EDA/BLYP-D3(BJ) MAE MSE MAX RMSE R2

DFT-SAPT

XB (88)d

Pni (6)

HB (23)

Disp (23)

Mix (20)

XB (88)

Pni (6)

HB (23)

Disp (23)

Mix (20)

1.02 −0.73 90 1.62 0.913

2.17 −2.17 61 2.33 0.927

0.38 −0.32 13 0.47 0.997

0.29 −0.29 21 0.33 0.993

0.16 −0.12 17 0.21 0.956

0.89 −0.46 88 2.14 0.887

0.92 −0.92 29 1.03 0.968

0.24 0.16 6 0.34 0.998

0.53 −0.52 35 0.58 0.984

0.27 −0.27 14 0.30 0.977

a

Mean absolute error (MAE, kcal/mol), mean signed error (MSE, kcal/mol), largest error relative to the interaction energy (MAX, %), root-meansquare error (RMSE, kcal/mol), coefficient of determination (R2). bIn the case of pnicogen-bonded complexes the reference data are calculated at the CCSD(T)/aug-cc-pVDZ level. cFor halogen-bonded (XB), pnicogen-bonded (Pni), hydrogen-bonded (HB), dispersion-dominated (Disp), and mixed (Mix) complexes. dNumber in parentheses corresponds to the number of complexes taken into consideration in a particular data set.

Table 2. Coefficients of Determination for Correlations between Particular Energy Terms of DFT-SAPT and EDA Methodsa XB Pni HB Disp Mix

repulsion

electrostatic

induction

dispersion

ind + disp

total

0.994 0.982 0.999 0.884 0.873

0.998 0.995 0.999 0.975 0.950

0.936 0.989 0.998 0.565 0.956

0.743 0.985 0.706 0.970 0.939

0.968 0.993 0.997 0.945 0.949

0.939 0.983 0.994 0.989 0.944

a

Rows represent different data sets: XB = halogen-bonded complexes; Pni = pnicogen-bonded complexes; HB = hydrogen-bonded complexes (from S66 data set); Disp = dispersion dominated complexes (from S66 data set); mix = “mixed” complexes (from S66 data set).

S66 Data Set. A set of complexes consists of 66 structures which can be divided into three types: hydrogen-bonded (HB, 23 complexes), dispersion-dominated (Disp, 23 complexes), and “other” (Mix, 20 complexes). 79 The DFT-SAPT interaction energy decomposition was performed with PBE0AC and aug-cc-pVDZ basis set. The dispersion term was scaled by a factor of 1.193.61 The interaction energy was recalculated accordingly to new values of the dispersion term.

used to estimate the CBS value from a single calculation. Namely, f = 1.193 is an average scaling factor for aug-cc-pVDZ calculation. The induction part summed with higher-order correction to (2) HF in DFT-SAPT, i.e., sum (E(2) ind + Eexch−ind + δ(HF)), contains information about intrafragment charge polarization and interfragment charge-transfer energy. Since separation of charge transfer from the remaining polarization is not an easy task, this term can be also consider as the upper bound to the so-called charge-transfer energy. This term will be compared with the corresponding orbital interaction energy term within the canonical EDA, designated as ΔEoi. Studied Data Sets. Halogen Bond (XB). A set of complexes contains 18 complexes from the X40 data set,62 46 complexes from the XB51 data set,63 8 complexes64 of crystal motifs taken from the Cambridge Structure Database, 17 structures of organic crystals taken from refs 65−70, and other complexes from refs 71−75. The structures of the halogen-bonded complexes were used from the original references without any additional optimization. DFT-SAPT calculations with density fitting2 were performed using PBE0/aug-cc-pVTZ (Tables 1− 6 in ref 76) and PBE0/aug-cc-pVDZ (Tables 7 and 8 in ref 76) methods. For large complexes (Tables 7 and 8), DFT-SAPT/ aug-cc-pVDZ dispersion energy results were scaled to aug-ccpVTZ basis set. The DFT-SAPT decomposition results were taken from ref 76. Pnicogen Bond (Pni). A set of complexes contains structures of benzene and its methyl and fluoro derivatives with halides of As, Sb,77 and Bi. The structures were optimized at DFT-D3/BLYP/def2-TZVPP level of theory with the Turbomole package. The DFT part of the DFT-SAPT calculations was treated using the LPBE0AC functional and aug-cc-pVDZ basis set. Then the dispersion energy was scaled to aug-cc-pVTZ basis set by a factor of 1.14,77 which is a ratio between the dispersion energies of the benzene···SbCl3 complex determined using aug-cc-pVDZ and aug-cc-pVTZ basis sets, and the interaction energy was recalculated.78

■

RESUTS AND DISCUSSION

Total Interaction Energies. First, we compared the performance of DFT-SAPT/PBE0AC/aug-cc-pVTZ and DFT-D3(BJ)/BLYP/TZ2P methods evaluating their accuracy in calculation of the total interaction energies for different types of noncovalent interactions. Table 1 shows statistical data for both applied methods with reference to CCSD(T)/CBS reference results, obtained by simple linear regression analysis. As can be seen, the applied methods yield well-correlated results, with coefficient of determination (R2) greater than 0.9 for all types of interactions, whereas their accuracy is higher for HB, Disp, and Mix complexes from the S66 data set and lower for XB and Pni complexes. The larger errors for halogenbonded complexes is associated with a broader data set which covers various interaction types according to their nature (electrostatic, charge-transfer, dispersion, or their combination) and wide range of stabilization energies. While the X40 data set provides accurate values for model complexes with RMSE = 0.39 kcal/mol at BLYP-D3(BJ)/def2-QZVP level of theory,62 the XB51 benchmark set yields worse results for GGA-DFT functionals with RMSE = 1.74 kcal/mol at the BLYP/aVTZ+CP level, and Grimme’s dispersion correction does not generally improve the performance.63 For pnicogenbonded complexes, both methods also give quite large RMSE values, but for such types of complexes only CCSD(T)/aug-ccpVDZ results for 6 complexes are available. Unfortunately, some specific classes of interactions such as the pnicogen bond are not included in the universal sets,80 and we cannot verify our finding. A small data set with the pnicogen bond in very D


Article

Journal of Chemical Theory and Computation small model systems proposed by Bauzá et al.81 gives MAE values of 0.96 kcal/mol for BLYP-D3 and 1.08 kcal/mol for PBE0 functionals with aug-cc-pVTZ basis set. Taking into account mean absolute error (MAE) and rootmean-square error (RMSE), which characterize the overall accuracy of methods, we can conclude that both methods provide quite comparable accuracy. This outcome allows us to compare different energetic terms obtained by two decomposition schemes, DFT-SAPT and EDA. Comparison of the Energy Components: DFT-SAPT vs Canonical EDA. In this part, the coefficients of determination between particular energy components will be discussed for all investigated data sets (Table 2). More statistical parameters for linear regression model EEDA = a· EDFT‑SAPT + b are collected in Table 2SI. The following energy terms of DFT-SAPT and canonical EDA schemes will be compared: Electrostatic: (1) Eelst vs ΔVelstat

Figure 1. Principal structures of halogen-bonded complexes used in this work (Hal = Cl, Br, or I).

Repulsion: (1) Eexch vs ΔEPauli

Induction: The difference between DFT-D3 dispersion correction and DFT-SAPT dispersion energy can be possibly attributed to the double counting of correlation effects in the middle-range region. To test this assumption we have calculated two more D3-based empirical dispersion corrections parametrized for the Hartree−Fock model51 as well as for the DFT-SAPT dispersion energy.17 The results showed a notable effect of double counting on the correlation between dispersion energy components. We have obtained R2 = 0.941 for the HF model and R2 = 0.955 for the DFT-SAPT model. This means that the aforementioned dispersion correction correlates significantly higher with the DFT-SAPT dispersion term than D3(BJ) correction (R2 = 0.743). Improved coefficients of determination are comparable with R2 for other energy components for a whole set of halogen bonded complexes. The repulsion energy term from EDA shows very good correlation with DFT-SAPT results but with systematically larger values for EDA than for DFT-SAPT. This difference is larger for complexes with short distance between interacting atoms, for which the repulsion term is very large. Since the main purpose of the energy decomposition is to gain insight into the nature of the interactions, we will next compare the relative contributions of individual energy terms to the total interaction energy. Among 122 halogen-bonded complexes, for 18 complexes two schemes determined a different nature of the interaction. According to canonical EDA scheme, the prevalent energy term in 8 complexes is the electrostatic interaction and in 10 complexes it is the induction term, whereas taking into account the DFT-SAPT method, the interactions in all aforementioned complexes are dispersiondominated. It should be emphasized that these discrepancies concern only very weakly interacting systems with total interaction energy weaker than −3 kcal/mol and with almost equal importance of different energy components. This is not a surprise, because D3 correction usually underestimates the total dispersion contribution and even minor difference in absolute values can change the conclusion about the nature of bonding in such complexes.

(2) (2) Eind = (Eind + Eexch − ind + δ(HF )) vs ΔEoi

Dispersion: (2) (2) Edisp = (Edisp + Eexch − disp) vs ΔEdisp ,

(Eind + Edisp)DFT − SAPT vs (ΔEoi + ΔEdisp)EDA

Total energy: DFT − SAPT EDA Etot vs Etot

Halogen Bond Data Set. The halogen-bonded complexes investigated in ref 76 can be split into two groups based on their stabilization energies: (1) strong complexes with a large charge-transfer contribution and stabilization energy in the range 7−32 kcal/mol and (2) complexes where dispersion energy is mostly dominant and with stabilization energies smaller than 7 kcal/mol. Figure 1 represents the main types of the studied halogen-bonded complexes, whereas Figure 2 shows the complexes with dihalogen bonds included also in the data set. The total interaction energies obtained by canonical EDA and DFT-SAPT methods are well correlated (R2 = 0.939). If we compare particular energy components (Figure 3), we can see that all components correlate well (R2 > 0.93), except for the dispersion term (R2 = 0.743). The outliers for induction and dispersion terms are two charge-transfer complexes, namely, diazabicyclo[2.2.2]octane (DABCO) or 1,3-dithiole2-thione-4-carboxylic acid (DTCA) interacting with diiodine (Figure 4). Such large differences between the D3 correction and the DFT-SAPT dispersion term were already reported for these complexes.82 In turn, underestimated dispersion in EDA is compensated by an overestimated induction term with respect to DFT-SAPT counterpart. At the same time, a sum of (ΔEoi + ΔEdisp)EDA and (Eind + Edisp)DFT−SAPT for all complexes are consistent with each other (R2 = 0.968). This means that in the EDA scheme a part of dispersion is included in the orbital interaction term. E


Article


Figure 2. Structures with dihalogen bonds: (a, b) o-diodoperfluorobenzene dimer; (c) 3,4,5-trichlorophenol dimer.

Figure 3. Correlations for total interaction energy and its different components (in kcal/mol) obtained by EDA and DFT-SAPT for complexes with a halogen bond.

complexes are dispersion dominated, and particularly in this case, the Grimme’s D3 correction totally coincides with the dispersion term obtained using the DFT-SAPT method (R2 = 0.985). However, as in the previous case, even better correlation can be observed between (ΔEoi + ΔEdisp)EDA and (Eind + Edisp)DFT−SAPT sums (R2 = 0.993). According to the nature of interactions, only 4 of 36 pnicogen-bonded complexes were found differently. All of them belong to the family of complexes with benzene. In complexes with benzene and hexamethylbenzene, the nature of the interaction is determined by the dispersion and electrostatic attraction almost equally (Table 3SI);77 therefore, assignment of prevalent significance to one or another energy term is not possible. In turn, in complexes with hexafluorobenzene, the dominant role of dispersion (almost 60% of all bonding forces) was uniquely determined by both methods. S66 Data Set. Hydrogen-Bonded Complexes. Twentythree hydrogen-bonded complexes in the S66 data set cover a range of interaction energy from −3 to −19 kcal/mol (Figure 7). The coefficient of determination for total energy values determined by the EDA and DFT-SAPT methods is the highest among all complexes studied (R2 = 0.994). The best correlations were also found for individual terms of energy. The only exception is the dispersion term with R2 = 0.706. For all H-bonded complexes, an electrostatic term was found to be the dominant one by both energy decomposition schemes. The electrostatic energy represents more than half of the total attraction force, i.e., sum of electrostatic, induction/ orbital interaction, and dispersion components. In turn, the dispersion estimated by the DFT-SAPT procedure is comparable in magnitude with the induction/orbital interaction term or is even greater. It should be noted that these values are always greater than D3 correction due to rather short H-bond distances and missing short-range dispersion in the Grimme’s correction.

Figure 4. Structures of complexes between diiodine and diazabicyclo[2.2.2]octane (DABCO) or 1,3-dithiole-2-thione-4-carboxylic acid (DTCA).

Pnicogen Bond. The principal structures of pnicogenbonded complexes investigated in this work are presented in Figure 5. There are two families of complexes: (i) with

Figure 5. Principal structures of complexes with pnicogen bond used in this work: (a) with benzene; (b) with hexamethylbenzene; and (c) with hexafluorobenzene. Pni = As, Sb, Bi; Hal = F, Cl, Br, I.

benzene and its hexamethyl derivative and (ii) with perfluorobenzene, which differ by various orientations of the PniHal3 molecule (Pni = As, Sb, Bi; Hal = F, Cl, Br, I). The total interaction energies as well as their appropriate components obtained by EDA and DFT-SAPT methods are perfectly correlated (R2 > 0.98) (Figure 6). Most likely it could be ascribed to similarity in their structures. Most of these F


Article


Figure 6. Correlations for total interaction energy and its different components (in kcal/mol) obtained by EDA and DFT-SAPT for complexes with pnicogen bond.

Figure 7. Correlations for total interaction energy and its different energy components (in kcal/mol) obtained by EDA and DFT-SAPT for Hbonded complexes from the S66 data set.

Figure 8. Correlations for total interaction energy and its different components (in kcal/mol) obtained by EDA and DFT-SAPT for dispersiondominated complexes from the S66 data set.

Dispersion-Dominated Complexes. This set consists of complexes with very weak or weak noncovalent interactions with energies from −1 to −10 kcal/mol. Total DFT-SAPT and EDA interaction energies correlate highly; the coefficient of determination is 0.989. Particular energy terms are also quite well correlated (Figure 8). The only exception is correlation between the induction part from DFT-SAPT and the orbital interaction term from canonical EDA with R2 = 0.565. Most likely such a low coefficient of determination can be explained

by small magnitude of these energy terms as expected for dispersion-dominated complexes and their narrow range (Table 2SI). For this same reason there is a more narrow range of data results in a slightly smaller coefficient of determination than in other types of complexes for the repulsion energy term. DFT-SAPT results for dispersion are regularly slightly larger than D3 correction and represent ca. 60−70% of total attractive forces. G


Article


Figure 9. Correlations for total interaction energy and its different components (in kcal/mol) obtained by EDA and DFT-SAPT for mixed complexes from the S66 data set.

Mixed Complexes. The last 20 complexes of the S66 data set contain other types of the most common noncovalent interactions in biomolecules (e.g., X-H···π (X = C, O, N) interactions, T-shaped aromatic ring complexes, nonspecific interactions of polar molecules) with various contributions of electrostatic, induction, and dispersion energetic terms. Taking into acccount their different nature, the coefficients of determination for total interaction energies and their particular terms are quite good. The smallest coefficient of determination was found for the repulsion term (R2 = 0.873). Its smaller magnitude compared to the other types of complexes can be explained by a more narrow range of energy data (Table 2SI). Nevetheless, this value is high enough to make a conclusion about good correlation. Both energy decomposition schemes provide very similar results (Figure 9). D3 dispersion correction is regularly smaller than the DFT-SAPT dispersion contribution. The nature of interactions is almost equally identified by both methods. Only systems where the contributions of the electrostatic component and the dispersion term are comparable in magnitude (e.g., benzene−water complex with OH···π contact) have small difficulty. Additionally, we have analyzed a compatibility of the decomposition schemes for the whole S66 data set (Figure 4SI). As in previous cases, the correlations between interaction energy as well as energy components are perfect (R2 > 0.98) with the exception of dispersion contribution (Figure 5SI). Detailed examination of correlation for the dispersion term shows that D3 correction values are always smaller than DFTSAPT dispersion. Moreover, the outliers are five complexes with a double hydrogen bond, i.e., uracil−uracil, AcOH− AcOH, AcNH2−AcNH2, AcOH−uracil, AcNH2−uracil, probably because of the even shorter distance between the interacting monomers.

increasingly important and a residual supermolecular Hartree− Fock term does not provide as good accuracy of the interaction energies as in conventional supermolecular methods. Second, correlations between corresponding energy terms of both methods are notably high, and in most cases both methods explain the character of noncovalent bonding in a very similar manner. Hence, both methods investigated in this work represent comparably valuable tools for analysis of noncovalent chemical bonds. In cases where the use of DFT-SAPT method is not feasible either due to computational complexity or due to limitations of the perturbation theory, we can recommend utilizing canonical EDA.

■

ASSOCIATED CONTENT

S Supporting Information *

The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.jctc.8b00034.

■

Testing of BLYP-D3 and PBE0-D3 functionals (Table 1SI, Figures 1SI−3SI); statistical data for energy term correlations (Table 2SI); percentage contribution of energetic terms for complexes with a pnicogen bond (Table 3SI); and correlations for interaction energy and its different components for S66 data set (Figures 4SI and 5SI) (PDF)

AUTHOR INFORMATION

Corresponding Author

*(P.H.) E-mail: [email protected]. ORCID

Robert Sedlak: 0000-0002-9900-5965 Pavel Hobza: 0000-0001-5292-6719

■

Notes

CONCLUSIONS The presented study shows that both decomposition schemes, DFT-SAPT and EDA, in combination with their most frequently used functionals and basis sets provide very similar results for the determination of the nature of noncovalent intermolecular interactions. First, their overall accuracy in the calculation of the total interaction energies is comparable which allows us to correlate particular energy components. Discrepancies appear only for complexes with short intermolecular distances, where higher-order effects become

The authors declare no competing financial interest.

■

ACKNOWLEDGMENTS The authors would like to thank Dr. Pecina for inspirational conversations which helped to shape this study in its early stage. O.A.S., R.S., and P.H. gratefully acknowledge support from the Czech Science Foundation (P208/12/G016). C.F.G. thanks The Netherlands Organization for Scientific Research (NWO/CW) for financial support. H


Article


■

(21) Fedorov, D. G.; Kitaura, K. Extending the Power of Quantum Chemistry to Large Systems with the Fragment Molecular Orbital Method. J. Phys. Chem. A 2007, 111, 6904−6914. (22) Ziegler, T.; Rauk, A. On the calculation of bonding energies by the Hartree Fock Slater method. Theor. Chim. Acta 1977, 46, 1−10. (23) Raupach, M.; Tonner, R. A periodic energy decomposition analysis method for the investigation of chemical bonding in extended systems. J. Chem. Phys. 2015, 142, 194105. (24) Langlet, J.; Bergès, J.; Reinhardt, P. Decomposition of intermolecular interactions: comparison between SAPT and densityfunctional decompositions. J. Mol. Struct.: THEOCHEM 2004, 685, 43−56. (25) Vyboishchikov, S. F.; Krapp, A.; Frenking, G. Two complementary molecular energy decomposition schemes: The Mayer and Ziegler−Rauk methods in comparison. J. Chem. Phys. 2008, 129, 144111. (26) Flick, J. C.; Kosenkov, D.; Hohenstein, E. G.; Sherrill, C. D.; Slipchenko, L. V. Accurate Prediction of Noncovalent Interaction Energies with the Effective Fragment Potential Method: Comparison of Energy Components to Symmetry-Adapted Perturbation Theory for the S22 Test Set. J. Chem. Theory Comput. 2012, 8, 2835−2843. (27) Phipps, M. J. S.; Fox, T.; Tautermann, C. S.; Skylaris, C.-K. Energy decomposition analysis approaches and their evaluation on prototypical protein-drug interaction patterns. Chem. Soc. Rev. 2015, 44, 3177−3211. (28) Schneider, W. B.; Bistoni, G.; Sparta, M.; Saitow, M.; Riplinger, C.; Auer, A. A.; Neese, F. Decomposition of Intermolecular Interaction Energies within the Local Pair Natural Orbital Coupled Cluster Framework. J. Chem. Theory Comput. 2016, 12, 4778−4792. (29) Kitaura, K.; Morokuma, K. A new energy decomposition scheme for molecular interactions within the Hartree-Fock approximation. Int. J. Quantum Chem. 1976, 10, 325−340. (30) Chen, W.; Gordon, M. S. Energy Decomposition Analyses for Many-Body Interaction and Applications to Water Complexes. J. Phys. Chem. 1996, 100, 14316−14328. (31) Glendening, E. D.; Streitwieser, A. Natural energy decomposition analysis: An energy partitioning procedure for molecular interactions with application to weak hydrogen bonding, strong ionic, and moderate donor-acceptor interactions. J. Chem. Phys. 1994, 100, 2900−2909. (32) Glendening, E. D. Natural Energy Decomposition Analysis: Explicit Evaluation of Electrostatic and Polarization Effects with Application to Aqueous Clusters of Alkali Metal Cations and Neutrals. J. Am. Chem. Soc. 1996, 118, 2473−2482. (33) Glendening, E. D. Natural Energy Decomposition Analysis: Extension to Density Functional Methods and Analysis of Cooperative Effects in Water Clusters. J. Phys. Chem. A 2005, 109, 11936−11940. (34) Werner, H.-J.; Knowles, P. J.; Knizia, G.; Manby, F. R.; Schütz, M.; Celani, P.; Korona, T.; Lindh, R.; Mitrushenkov, A.; Rauhut, G.; Shamasundar, K. R.; Adler, T. B.; Amos, R. D.; Bernhardsson, A.; Berning, A.; Cooper, D. L.; Deegan, M. J. O.; Dobbyn, A. J.; Eckert, F.; Goll, E.; Hampel, C.; Hesselmann, A.; Hetzer, G.; Hrenar, T.; Jansen, G.; Köppl, C., Liu, Y.; Lloyd, A. W.; Mata, R. A.; May, A. J.; McNicholas, S. J.; Meyer, W.; Mura, M. E.; Nicklass, A.; O’Neill, D. P.; Palmieri, P.; Pflüger, K.; Pitzer, R.; Reiher, M.; Shiozaki, T.; Stoll, H.; Stone, A. J.; Tarroni, R.; Thorsteinsson, T.; Wang, M.; Wolf, A. MOLPRO, version 2010.1, a package of ab initio programs. http:// www.molpro.net. (35) SCM. ADF 2014; Theoretical Chemistry, Vrije Universiteit: Amsterdam, The Netherlands. http://www.scm.com. (36) Heßelmann, A.; Jansen, G. First-order intermolecular interaction energies from Kohn−Sham orbitals. Chem. Phys. Lett. 2002, 357, 464−470. (37) Heßelmann, A.; Jansen, G. Intermolecular induction and exchange-induction energies from coupled-perturbed Kohn−Sham density functional theory. Chem. Phys. Lett. 2002, 362, 319−325.

REFERENCES

(1) Wu, Q.; Ayers, P. W.; Zhang, Y. Density-based energy decomposition analysis for intermolecular interactions with variationally determined intermediate state energies. J. Chem. Phys. 2009, 131, 164112. (2) Hesselmann, A.; Jansen, G.; Schutz, M. Density-functional theory-symmetry-adapted intermolecular perturbation theory with density fitting: A new efficient method to study intermolecular interaction energies. J. Chem. Phys. 2005, 122, 014103. (3) Fedorov, D. G.; Kitaura, K. Pair interaction energy decomposition analysis. J. Comput. Chem. 2007, 28, 222−237. (4) Su, P.; Li, H. Energy decomposition analysis of covalent bonds and intermolecular interactions. J. Chem. Phys. 2009, 131, 014102. (5) Horn, P. R.; Mao, Y.; Head-Gordon, M. Probing non-covalent interactions with a second generation energy decomposition analysis using absolutely localized molecular orbitals. Phys. Chem. Chem. Phys. 2016, 18, 23067−23079. (6) Bickelhaupt, F. M.; Baerends, E. J. Kohn-Sham Density Functional Theory: Predicting and Understanding Chemistry. In Reviews in Computational Chemistry; Lipkowitz, K. B., Boyd, D. B., Eds.; Wiley: New York, 2000; Vol. 15, pp 1−86. (7) Morokuma, K. Molecular Orbital Studies of Hydrogen Bonds. III. C = O···H−O Hydrogen Bond in H2CO···H2O and H2CO··· 2H2O. J. Chem. Phys. 1971, 55, 1236. (8) Jeziorski, B.; Moszynski, R.; Szalewicz, K. Perturbation Theory Approach to Intermolecular Potential Energy Surfaces of van der Waals Complexes. Chem. Rev. 1994, 94, 1887−1930. (9) Khaliullin, R. Z.; Cobar, E. A.; Lochan, R. C.; Bell, A. T.; HeadGordon, M. Unravelling the origin of intermolecular interactions using absolutely localized molecular orbitals. J. Phys. Chem. A 2007, 111, 8753−8765. (10) Glendening, E. D.; Streitwieser, A. Natural energy decomposition analysis: An energy partitioning procedure for molecular interactions with application to weak hydrogen bonding, strong ionic, and moderate donor−acceptor interactions. J. Chem. Phys. 1994, 100, 2900−2909. (11) Steinmann, S. N.; Corminboeuf, C.; Wu, W.; Mo, Y. Dispersion-Corrected Energy Decomposition Analysis for Intermolecular Interactions Based on the BLW and dDXDM Methods. J. Phys. Chem. A 2011, 115, 5467−5477. (12) Fonseca Guerra, C.; Bickelhaupt, F. M.; Snijders, J. G.; Baerends, E. J. The Nature of the Hydrogen Bond in DNA Base Pairs: The Role of Charge Transfer and Resonance Assistance. Chem. - Eur. J. 1999, 5, 3581−3594. (13) Stone, A. J.; Misquitta, A. J. Charge-transfer in symmetryadapted perturbation theory. Chem. Phys. Lett. 2009, 473, 201−205. (14) Podeszwa, R. Interactions of graphene sheets deduced from properties of polycyclic aromatic hydrocarbons. J. Chem. Phys. 2010, 132, 044704. (15) Hesselmann, A.; Korona, T. Intermolecular symmetry-adapted perturbation theory study of large organic complexes. J. Chem. Phys. 2014, 141, 094107. (16) Hesselmann, A. Comparison of Intermolecular Interaction Energies from SAPT and DFT Including Empirical Dispersion Contributions. J. Phys. Chem. A 2011, 115, 11321−11330. (17) Sedlak, R.; Ř ezać, J. Empirical D3 Dispersion as a Replacement for ab Initio Dispersion Terms in Density Functional Theory-Based Symmetry-Adapted Perturbation Theory. J. Chem. Theory Comput. 2017, 13, 1638−1646. (18) Church, J.; Pezeshki, S.; Davis, C.; Lin, H. Charge Transfer and Polarization for Chloride Ions Bound in CIC Transport Proteins: Natural Bond Orbital and Energy Decomposition Analyses. J. Phys. Chem. B 2013, 117, 16029−16043. (19) Thellamurege, N.; Hirao, H. Water Complexes of Cytochrome P450: Insights from Energy Decomposition Analysis. Molecules 2013, 18, 6782−6791. (20) Hirao, H. Energy decomposition analysis of the protein environmental effect: The case of cytochrome P450cam compound I. Chem. Lett. 2011, 40, 1179−1181. I


Article

Journal of Chemical Theory and Computation (38) Heßelmann, A.; Jansen, G. Intermolecular dispersion energies from time-dependent density functional theory. Chem. Phys. Lett. 2003, 367, 778−784. (39) Ziegler, T.; Rauk, A. Carbon monoxide, carbon monosulfide, molecular nitrogen, phosphorus trifluoride, and methyl isocyanide as σ donors and π acceptors. A theoretical study by the Hartree-FockSlater transition-state method. Inorg. Chem. 1979, 18, 1755−1759. (40) Ziegler, T.; Rauk, A. A theoretical study of the ethylene-metal bond in complexes between Cu+, Ag+, Au+, Pt0 or Pt2+ and ethylene, based on the Hartree-Fock-Slater transition-state method. Inorg. Chem. 1979, 18, 1558−1565. (41) Jansen, G. Symmetry-adapted perturbation theory based on density functional theory for noncovalent interactions. WIREs Comput. Mol. Sci. 2014, 4, 127−144. (42) te Velde, G.; Bickelhaupt, F. M.; Baerends, E. J.; Fonseca Guerra, C.; van Gisbergen, S. J. A.; Snijders, J. G.; Ziegler, T. Chemistry with ADF. J. Comput. Chem. 2001, 22, 931−967. (43) Misquitta, A. J.; Podeszwa, R.; Jeziorski, B.; Szalewicz, K. Intermolecular potentials based on symmetry-adapted perturbation theory including dispersion energies from time-dependent density functional calculations. J. Chem. Phys. 2005, 123, 214103. (44) Hohenstein, E. G.; Sherrill, C. D. Density fitting of intramonomer correlation effects in symmetry-adapted perturbation theory. J. Chem. Phys. 2010, 133, 014101. (45) Korona, T. A coupled cluster treatment of intramonomer electron correlation within symmetry-adapted perturbation theory: benchmark calculations and a comparison with a density-functional theory description. Mol. Phys. 2013, 111, 3705−3715. (46) Grüning, M.; Gritsenko, O. V.; van Gisbergen, S. J. A.; Baerends, E. J. Shape corrections to exchange-correlation potentials by gradient-regulated seamless connection of model potentials for inner and outer region. J. Chem. Phys. 2001, 114, 652. (47) TURBOMOLE, V6.3 2011, a development of the University of Karlsruhe and the Forschungszentrum Karlsruhe GmbH, 1989−2007, TURBOMOLE GmbH, since 2007. http://www.turbomole.com. (48) Bickelhaupt, F. M.; Nibbering, N. M. M.; van Wezenbeek, E. M.; Baerends, E. J. Central bond in the three CN· dimers NC-CN, CN-CN and CN-NC: electron pair bonding and Pauli repulsion effects. J. Phys. Chem. 1992, 96, 4864−4873. (49) Bickelhaupt, F. M.; Diefenbach, A.; de Visser, S. P.; de Koning, L. J.; Nibbering, N. M. M. Nature of the Three-Electron Bond in H2SSH2+. J. Phys. Chem. A 1998, 102, 9549−9553. (50) Grimme, S.; Antony, J.; Ehrlich, S.; Krieg, H. A consistent and accurate ab initio parametrization of density functional dispersion correction (DFT-D) for the 94 elements H-Pu. J. Chem. Phys. 2010, 132, 154104. (51) Grimme, S.; Ehrlich, S.; Goerigk, L. Effect of the damping function in dispersion corrected density functional theory. J. Comput. Chem. 2011, 32, 1456−1465. (52) Becke, A. D. Density-functional exchange-energy approximation with correct asymptotic-behavior. Phys. Rev. A: At., Mol., Opt. Phys. 1988, 38, 3098−3100. (53) Lee, C.; Yang, W.; Parr, R. G. Development of the ColleSalvetti correlation-energy formula into a functional of the electrondensity. Phys. Rev. B: Condens. Matter Mater. Phys. 1988, 37, 785−789. (54) Goerigk, L.; Grimme, S. A thorough benchmark of density functional methods for general main group thermochemistry, kinetics, and noncovalent interactions. Phys. Chem. Chem. Phys. 2011, 13, 6670−6688. (55) Fonseca Guerra, C.; van der Wijst, T.; Swart, M.; Poater, J.; Bickelhaupt, F. M. Adenine versus Guanine Quartets in Aqueous Solution. Dispersion-Corrected DFT Study on the Differences in πStacking and Hydrogen-Bonding Behaviour. Theor. Chem. Acc. 2010, 125, 245−252. (56) Brauer, B.; Kesharwani, M. K.; Kozuch, S.; Martin, J. M. L. The S66 × 8 benchmark for noncovalent interactions revisited: explicitly correlated ab initio methods and density functional theory. Phys. Chem. Chem. Phys. 2016, 18, 20905−20925.

(57) DiLabio, G. A.; Otero-de-la-Roza, A. Noncovalent Interactions in Density Functional Theory. In Reviews in Computational Chemistry; Parrill, A. L., Lipkowitz, K. B., Eds.; John Wiley & Sons: Hoboken, NJ, 2016; Vol. 29, pp 1−97. (58) Snijders, J. G.; Baerends, E. J.; Vernooijs, P. Roothaan-HartreeFock-Slater atomic wave functions. Single-zeta, double-zeta, and extended Slater-type basis sets for 87Fr-103Lr. At. Data Nucl. Data Tables 1981, 26, 483−509. (59) van Lenthe, E.; Baerends, E. J.; Snijders, J. G. Relativistic regular two-component Hamiltonians. J. Chem. Phys. 1993, 99, 4597− 4610. (60) van Lenthe, E.; Baerends, E. J.; Snijders, J. G. Relativistic total energy using regular approximations. J. Chem. Phys. 1994, 101, 9783− 9792. (61) Rezac, J.; Hobza, P. Extrapolation and Scaling of the DFTSAPT Interaction Energies toward the Basis Set Limit. J. Chem. Theory Comput. 2011, 7, 685−689. (62) Rezac, J.; Riley, K. E.; Hobza, P. Benchmark Calculations of Noncovalent Interactions of Halogenated Molecules. J. Chem. Theory Comput. 2012, 8, 4285−4292. (63) Kozuch, S.; Martin, J. M. L. Halogen Bonds: Benchmarks and Theoretical Analysis. J. Chem. Theory Comput. 2013, 9, 1918−1931. (64) Wang, W. Halogen Bond Involving Hypervalent Halogen: CSD Search and Theoretical Study. J. Phys. Chem. A 2011, 115, 9294− 9299. (65) Trnka, J.; Sedlak, R.; Kolar, M.; Hobza, P. Differences in the Sublimation Energy of Benzene and Hexahalogenbenzenes Are Caused by Dispersion Energy. J. Phys. Chem. A 2013, 117, 4331− 4337. (66) Mukherjee, A.; Desiraju, G. R. Halogen Bonding and Structural Modularity in 2,3,4- and 3,4,5-Trichlorophenol. Cryst. Growth Des. 2011, 11, 3735−3739. (67) Cincic, D.; Friscic, T.; Jones, W. A Stepwise Mechanism for the Mechanochemical Synthesis of Halogen-Bonded Cocrystal Architectures. J. Am. Chem. Soc. 2008, 130, 7524−7525. (68) Cauliez, P.; Polo, V.; Roisnel, T.; Llusar, R.; Fourmigue, M. The thiocyanate anion as a polydentate halogen bond acceptor. CrystEngComm 2010, 12, 558−566. (69) Cincic, D.; Friscic, T.; Jones, W. Experimental and database studies of three-centered halogen bonds with bifurcated acceptors present in molecular crystals, cocrystals and salts. CrystEngComm 2011, 13, 3224−3231. (70) Jay, J. I.; Padgett, C. W.; Walsh, R. D. B.; Hanks, T. W.; Pennington, W. T. Noncovalent Interactions in 2-Mercapto-1methylimidazole Complexes with Organic Iodides. Cryst. Growth Des. 2001, 1, 501−507. (71) Munusamy, E.; Sedlak, R.; Hobza, P. On the nature of the stabilization of benzene···dihalogen and benzene···dinitrogen complexes: CCSD(T)/CBS and DFT-SAPT calculations. ChemPhysChem 2011, 12, 3253−3261. (72) Riley, K. E.; Rezac, J.; Hobza, P. Competition between halogen, dihalogen and hydrogen bonds in bromo- and iodomethanol dimers. J. Mol. Model. 2013, 19, 2879−2883. (73) Sedlak, R.; Deepa, P.; Hobza, P. Why Is the L-Shaped Structure of X2···X2 (X = F, Cl, Br, I) Complexes More Stable Than Other Structures? J. Phys. Chem. A 2014, 118, 3846−3855. (74) Lu, Y. X.; Zou, J. V.; Fan, J. C.; Zhao, W. N.; Jiang, Y. J.; Yu, Q. S. Ab initio calculations on halogen-bonded complexes and comparison with density functional methods. J. Comput. Chem. 2009, 30, 725−732. (75) Chudzinski, M. G.; Taylor, M. S. Correlations between Computation and Experimental Thermodynamics of Halogen Bonding. J. Org. Chem. 2012, 77, 3483−3491. (76) Kolár,̌ M. H.; Deepa, P.; Ajani, H.; Pecina, A.; Hobza, P. Characteristics of a σ-Hole and the Nature of a Halogen Bond. Top. Curr. Chem. 2014, 359, 1−25. (77) Lo, R.; Š vec, P.; Růzǐ čková, Z.; Růzǐ čka, A.; Hobza, P. On the nature of the stabilisation of the E···π pnicogen bond in the SbCl3··· toluene complex. Chem. Commun. 2016, 52, 3500−3503. J


Article

Journal of Chemical Theory and Computation (78) Lo, R. In preparation. (79) Rezac, J.; Riley, K. E.; Hobza, P. S66: A Well-balanced Database of Benchmark Interaction Energies Relevant to Biomolecular Structures. J. Chem. Theory Comput. 2011, 7, 2427−2438. (80) Ř ezać, J.; Hobza, P. Benchmark Calculations of Interaction Energies in Noncovalent Complexes and Their Applications. Chem. Rev. 2016, 116, 5038−5071. (81) Bauzá, A.; Alkorta, I.; Frontera, A.; Elguero, J. On the Reliability of Pure and Hybrid DFT Methods for the Evaluation of Halogen, Chalcogen, and Pnicogen Bonds Involving Anionic and Neutral Electron Donors. J. Chem. Theory Comput. 2013, 9, 5201− 5210. (82) Deepa, P.; Sedlak, R.; Hobza, P. On the origin of the substantial stabilisation of the electron-donor 1,3-dithiole-2-thione-4-carboxyclic acid···I2 and DABCO···I2 complexes. Phys. Chem. Chem. Phys. 2014, 16, 6679−6686.

K


Comparison of the DFT-SAPT and Canonical EDA Schemes for the

Recommend Documents