Revised Damping Parameters for the D3 ... - ACS Publications

May 20, 2016 - Center for Computational Molecular Science and Technology, School of Chemistry and Biochemistry and School of Computational. Science an...
2 downloads 7 Views 1MB Size
Letter pubs.acs.org/JPCL

Revised Damping Parameters for the D3 Dispersion Correction to Density Functional Theory Daniel G. A. Smith,*,† Lori A. Burns,*,‡ Konrad Patkowski,*,† and C. David Sherrill*,¶ †

Department of Chemistry and Biochemistry, Auburn University, Auburn, Alabama 36849, United States Center for Computational Molecular Science and Technology, School of Chemistry and Biochemistry, Georgia Institute of Technology, Atlanta, Georgia 30332-0400, United States ¶ Center for Computational Molecular Science and Technology, School of Chemistry and Biochemistry and School of Computational Science and Engineering, Georgia Institute of Technology, Atlanta, Georgia 30332-0400, United States ‡

S Supporting Information *

ABSTRACT: Since the original fitting of Grimme’s DFT-D3 damping parameters, the number and quality of benchmark interaction energies has increased significantly. Here, conventional benchmark sets, which focus on minimum-orientation radial curves at the expense of angular diversity, are augmented by new databases such as side chain−side chain interactions (SSI), which are composed of interactions gleaned from crystal data and contain no such minima-focused bias. Moreover, some existing databases such as S22×5 are extended to shorter intermolecular separations. This improved DFT-D3 training set provides a balanced description of distances, covers the entire range of interaction types, and at 1526 data points is far larger than the original training set of 130. The results are validated against a new collection of 6773 data points and demonstrate that the effect of refitting the damping parameters ranges from no change in accuracy (LC-ωPBE-D3) to an almost 2-fold decrease in average error (PBE-D3).

D

While density functionals that include nonlocal dispersion either in an atom-pairwise form16−19 or through a nonlocal correlation functional20,21perform generally quite well at long-range (the leading asymptotic constants C6 are usually accurate to better than 10%), at distances slightly shorter than the van der Waals minimum the accuracy of dispersioncorrected DFT deteriorates quickly. Large errors of DFT-D222 at separations slightly shorter than equilibrium (0.9·Req) have been observed already in ref 2; however, these issues went relatively unnoticed by the general community. Recently, two of us4 performed extensive benchmark calculations of complexes involving carbon dioxide and coronene-sized models of graphene and carbon nanotubes at distances down to 0.8·Req, observing a dramatic drop in DFT-D accuracy at short-range. Reference 4 attributed the short-range issues to two sources: the overestimation of exchange due to incorrect asymptotic behavior of the exchange-correlation potential and the inadequacy of damping functions commonly used in atompairwise dispersion expressions. As evidenced by the long-range corrected functional LC-ωPBE-D323 performing only somewhat better than standard functionals, the second issue is likely more important. It should be noted that the errors of DFT-D based approaches become quite severe already at distances around 0.8−0.9·Req, where the interaction energies are either

ispersion-including density functional theory (DFT) methods have become very common for the evaluation of intermolecular interaction energies and potential energy surfaces. Numerous studies have shown that these methods can obtain impressive accuracy at the van der Waals minima, with average errors down to about 0.2 kcal/mol (0.8 kJ/mol) in some cases.1 However, the accuracy near the minima is not always representative of the accuracy across the entire potential energy surface.1−5 Databases of accurate noncovalent interactions (NCI) have been invaluable for assessing the performance of both established approaches and new approximations. The early databases6−8 focused on the van der Waals minimum geometries. It was realized shortly thereafter5 that good performance of a DFT-based approach at the minimum does not necessarily translate to a similar performance at other intermolecular separations; for example, semilocal functionals, even when parametrized for noncovalent interactions,9,10 cannot reproduce the correct long-range behavior of interaction energies. 11 Therefore, many newer databases 2,12,13 are composed of radial curves (i.e., potential energies as a function of intermolecular separation R) passing through the van der Waals minimum. Only a few benchmark studies14,15 have included nonminimum angular configurations of varied interacting monomers, but there is a multitude of fulldimensional PES (potential energy surface) data for individual small complexes to assess the accuracy of DFT-based approaches at nonminimum orientations. © 2016 American Chemical Society

Received: April 12, 2016 Accepted: May 20, 2016 Published: May 20, 2016 2197

DOI: 10.1021/acs.jpclett.6b00780 J. Phys. Chem. Lett. 2016, 7, 2197−2203

Letter

The Journal of Physical Chemistry Letters Table 1. Datasets Utilized in the Training and Validation Setsa database

points

curves

Training CO2·PAH HBC6 NBC10extc,d S22×7c,d SSI500 X31×10e Validation ACHC BBI C2H4·NT CH4·PAHc CO2·NPHAC S66×10c SSI Water2510 Thermochemistry CONFd

1526 249 118 195 154 500 310 6773 54 100 75 405 96 660 f 2873 2510 52 52

114 45 6 10 22 b

31 148 6 b

15 45 16 66 b b

b

largestg

reference

description

27 6 12 19 20 18

4 12, 28 1, 28 2, 8, here unpub. 29

CO2 with PAHs the size of benzene through coronene and curved coronene dissociation curves of doubly hydrogen-bonded complexes dissociation curves of dispersion-bound complexes dissociation curves for a mix of hydrogen bonded and dispersion bonded complexes a representative subset of 500 structures from SSI dissociation curves of organic halide, halohydride, and halogen complexes

19 10 26 25 27 16 21 2

30 unpub. unpub. 3, 31, here 32 13, here unpub. 33, 34, here

rise, twist, slide, shift, roll, and tilt of adenine:cytosine nucleobase step peptide backbone−backbone complexes ethene with curved coronene methane with PAHs the size of benzene through coronene and curved coronene CO2 with nitrogen-doped polyheterocyclic aromatic compounds dissociation curves for a balanced mix of biomolecule NCI bonding motifs peptide side chain-side chain complexes water dimer PES

20

35−39, here

ACONF, CYCONF, PCONF-A, and SCONF-A conformational data sets

All benchmark datasets are of MP2/CBS + ΔCCSD(T)/aDZ quality or better. Details on the benchmark level of theory are given in the Supporting Information. bDatabase does not contain curves. cDatabase was extended to shorter ranges. dDatabase was recomputed at a higher level of theory, see Supporting Information, Table SIV. eThe X40x10 database with complexes containing iodine removed. fSSI contains 3380 complexes. The stated figure is lowered by 500 from the SSI500 fitting subset and 7 for which GGA functionals do not reliably converge. gThe largest number of heavy atoms in the data set. a

Figure 1. Top Panel: Ternary diagram comparison between the original -D3 fitting set and the current fitting set. Bottom Panel: Ternary diagram breakdown of the three major categories of points included in the validation set.

equilibrium configurations. Therefore, we are now in a much better position to revisit the parameter selection in the original DFT-D316 and DFT-D3(BJ)27 approaches. The ultimate DFT method for NCI should reproduce CCSD(T)-quality results at both short and long intermolecular distances for a very diverse set of interactions. Therefore, we require (i) all benchmark databases to be computed at the MP2/CBS + ΔCCSD(T)/aDZ (≈5% accuracy) level of theory or better, (ii) all dissociation curves to contain at least one positive point on the repulsive wall and to extend to distances longer than the equilibrium position, and (iii) all databases taken together to represent a diverse set of NCI as demonstrated by energy component analysis. Upon consideration of databases for training and validation roles, it was

still negative or only slightly positive (up to, say, 10 kcal/mol). Such mildly repulsive configurations are probed in molecular dynamics simulations at standard conditions of temperature and pressure and the errors of DFT at this range are likely to adversely influence the resulting quantities such as spectra and virial coefficients. Moreover, vastly different accuracy at different R leads to highly inaccurate gradients of the intermolecular potential (forces).24 The availability of multiple benchmark NCI databases including radial curves, along with the increase in computer power, has made it possible to evaluate and optimize DFTbased approaches on a much larger scale, up to exploring complete DFT parameter spaces.25,26 Moreover, the very recent development of the SSI benchmark database by two of us and collaborators extends the available benchmark set far beyond 2198

DOI: 10.1021/acs.jpclett.6b00780 J. Phys. Chem. Lett. 2016, 7, 2197−2203

Letter

The Journal of Physical Chemistry Letters

Figure 2. Top Panel: All functionals utilizing the original damping parameters. Bottom Panel: All functionals utilizing the new damping parameters (as signified by the letter M (modified) in the functional name). The light gray outlines represent the MCURE of the original damping parameters. The right-hand panels give the spread of MCURE values for all functionals.

contributions from the three attractive SAPT components: electrostatics, induction, and dispersion. The “minima cross sections” ternary plot demonstrates that these types of data sets primarily explore intermolecular interactions that are dispersion- and electrostatically dominated. While this region can be argued to be the most important, it can be clearly seen that the SSI and Water2510 sets explore a more diverse set of interaction archetypes. The lack of structures bound by pure induction (top right corner of a ternary diagram) is considered quite satisfactory, as only unusual systems such as H+−LiH could be found in this region. The original -D3 training set of size 130 contains 72 intermolecular interactions and 58 thermochemistry-based data points. Of those intermolecular interactions (for which SAPT can be computed), 66 are the S22 and S22+27 geometries (S22×5 at z = 1.0, 1.2, 1.5) plotted in Figure 1 (top left) in comparison to the current training set (top right). S22+ lacks any points shorter than equilibrium (a region where damping plays a substantial role), and its original benchmark does not meet the standard of accuracy adhered to in this work. Thus, the new training set is a significant improvement in quality, diversity, and range. As thermochemistry is not the focus of the current work, the CONF data set (52 of the original 58 thermochemical data points) is used only to ensure that the fitting procedure does not radically skew this type of data. Two fitting forms are utilized. The Becke−Johnson (BJ) damping27,43 with three parameters s8, α1, and α2 is

discovered that many established databases themselves do not meet these requirements. Benchmark dissociation curves are typically created by first finding the global or local minimum of a dimer. Then, more points are created by changing the distance between the R monomers by a fraction of the reduced distance z = R . These eq

types of databases are denoted “minima cross section” and comprise all databases employed in this study (Table 1) except for ACHC (partially), BBI, SSI, and Water2510. Selection of the correct set of z distances is of utmost importance to accurately capture the shape of the potential energy curve. For example, the S22×5 database was originally selected z as 0.9, 1.0, 1.2, 1.5, and 2.0, which sufficiently covers long-range (z ≥ 1) regimes. However, the short range (z < 1) is inadequate: for z = 0.9, only one curve has a positive point on the repulsive wall. By extending z to 0.8, 11 curves reach positive points, while by 0.7, all 22 do. This problem is not isolated to the S22×5 database: in total, the S22×5 (S22×7), S66×8 (S66×10), NBC10 (NBC10ext), and CH4·PAH data sets were extended to shorter ranges for this Letter; new names are denoted in parentheses. To verify that all types of interactions are represented, symmetry-adapted perturbation theory (SAPT)40,41 energy decomposition was computed for all data sets. The results are plotted in Figure 1 using “ternary” diagrams42 that represent the nature of an intermolecular interaction by relative 2199

DOI: 10.1021/acs.jpclett.6b00780 J. Phys. Chem. Lett. 2016, 7, 2197−2203

Letter

The Journal of Physical Chemistry Letters BJ Edisp =−

1 2

∑ ∑

sn

A ≠ B n = 6,8

CnAB n rAB + (α1·R 0AB + α2)n

range inaccuracies in DFT+D can be attributed to KS exchange being too repulsive.45−47 Figure 2 shows that, after refitting, the overall MCURE for every functional either improves or, for LC-ωPBE-D3M, does not change. The long-range accuracies are now more consistent, with errors between 4 and 12%. As the C6 coefficients used in the dispersion term are only accurate to about 8% on the average,16 the overall performance at longrange is quite satisfactory. The accuracy for the SSI data set is also increased significantly. Before refitting, errors ranged from 7 to 26%, and now they only range from 4 to 14%. The reduction in short-range MCURE is not always guaranteed; for the LC-ωPBE-D3 and BP86-D3(BJ) functionals, the errors become slightly worse. However, for the B2PLYP, B3LYP, and BLYP series, short-range MCURE’s significantly improve, from 14−28% before refitting, down to 8−14% after refitting. The most notable single functional improvement comes from the refitting of the popular PBE-D3 method, for which the overall MCURE changed from 19% to 10% and the MCURE for the SSI validation set decreased from 26% to 9%. For the BLYP and PBE family of functionals, the variability between the GGA, hybrid, and double-hybrid (or long-range corrected) forms has been reduced so that the overall MCURE is within 5% for a given functional family. This reduction in variability between functional forms allows for functional selection based on criterion other than intermolecular accuracy, such as computational expense. The quantities shown in Figure 2 represent an overview of the statistics. To obtain a better understanding of the improvement per data set, the MCURE for each data set is shown in Figure 3. The BLYP-D3 functional was chosen as an

(1)

The “zero damping” or Chai−Head-Gordon (CHG) form,16,44 CHG Edisp =−

1 2

∑ ∑ A ≠ B n = 6,8

sn

CnAB 1 n rAB 1 + 6(rAB/(sr , nR 0AB) + R 0ABβ)−αn (2)

has three parameters s8, sr,6, and β (the sr,8 values are held constant at 1, and the s6 parameter in eqs 1 and 2 is fixed at 0.64 for B2PLYP and 1 otherwise). The α6 and α8 parameters are set to 14 and 16, respectively, and not adjusted. The β parameter is introduced to CHG damping in this paper to give the same number of fitting parameters as BJ damping. The training and validation databases cover a large range of interaction energies, making a metric like mean unsigned error heavily favor complexes with large interaction energies. A statistical quantity like the mean unsigned relative error would be a better metric to measure the overall accuracy; however, this methodology runs into singularities whenever the PES crosses the zero line (which happens at short range and may also happen at long range). To circumvent this issue, a capped relative error (CRE) is introduced where the weight is capped at a certain value. Unfortunately, a single value as in ref 4 is not appropriate for the entire database. Instead, the cap needs to have a functional form to take into account the large range of interaction energies. ⎛E − E ⎞ ref ⎟ CRE = ⎜⎜ ⎟ · 100%, ⎝ Eweight ⎠

⎧ ξ|Eref − eq | ⎫ ⎬ Eweight = max⎨|Eref | , z3 ⎭ ⎩ (3)

ξ is a flexible dimensionless parameter that determines the severity of the capping. In this paper, a value of 0.2 was selected as it represents a good balance between the number of points capped and avoiding singularity. At this weight, only 10% of the values had their weights altered, with an average difference between Eweight and the benchmark (Eref) of 0.25 kcal/mol. For the SSI, BBI, and CONF data sets, complete curves are not available, and the z value is undefined. In this case, a simple cap of 0.5 kcal/mol was used; this value was chosen because the capping affects a similar percentage of overall points (7%). Variations of the CRE include capped unsigned relative error (CURE) and mean capped unsigned relative error (MCURE). For the overall statistical quantities, each database is weighted equally, except for SSI, which always contributes 1 to the 3 overall statistic. This was done so that the different number of points in each data set does not implicitly weight the data sets and the relative importance of SSI is always represented. All damping parameters were fitted to the training set and all statistical quantities are given on the validation set. Additional statistics can be found in the SI. The MCURE values for all functionals with original damping parameters are shown in the top panel of Figure 2. The diversity in the statistics between functionals and damping forms is quite striking, especially for the SSI and short-range quantities. For this large data set, it is apparent that long-range DFT+D is overall quite excellent, often with below 10% MCURE and at worst case 18% MCURE for the B97-D3(BJ) functional. However, the accuracy radically deteriorates at short-range, with a best case scenario of 12% for LC-ωPBE-D3 and worst case of 29% for B97-D3. Indeed, some of the short-

Figure 3. Performance of the BLYP-D3 (outlines) and BLYP-D3M (filled bars) functionals on separate subsets of the validation data set.

example, as it represents a reasonable middle ground in terms of overall improvement. Figure 3 demonstrates that the accuracy per data set can vary dramatically, especially at short range. After the refitting, we can observe that the accuracy per data set is much more even, and the largest outliers have been brought more in line with the remaining data sets. To examine the fundamental limits of the current damping forms, Figure 4 shows the CRE per data point. For CHG damping, dispersion-dominated systems are typically underbound while electrostatically bound systems are overbound. Interestingly, the refitting of the CHG damping appears to fix the electrostatically bound systems, while the errors remain relatively unchanged for dispersion-bound systems. The overlapping overbound and underbound points on the ternary 2200

DOI: 10.1021/acs.jpclett.6b00780 J. Phys. Chem. Lett. 2016, 7, 2197−2203

Letter

The Journal of Physical Chemistry Letters

intermolecular separation and the interaction type. Using this new data set, we revisited the determination of the DFT-D3 damping parameters, as their original values16,27 were not fitted to short range points and consequently do poorly in this region. Refitting can reduce DFT-D3 errors by up to half and in general greatly reduces the large variability between different functional selections. We recommend the B2PLYP, B3LYP, and BLYP family of functionals paired with the -D3M(BJ) dispersion correction due to their consistent accuracy at all ranges.



COMPUTATIONAL METHODS All computations were performed using the MOLPRO48 2012 and PSI449 suites of ab initio quantum chemistry programs. The B2PLYP,50 B3LYP,51,52 BLYP,53,54 BP86,53,55 PBE0,56,57 LCωPBE,23 PBE,58 and B9722 functionals were utilized. The def2QZVP59 and aug-cc-pVDZ60,61 (shortened to QZVP and aDZ, respectively) basis sets were utilized. Unless otherwise stated, all computations were counterpoise corrected62,63 and utilized the QZVP basis set. Density fitting was utilized for all DFT computations.64,65 The SAPT calculations were performed at the SAPT0/jun-cc-pVDZ66,67 or higher level using the PSI4 code.

Figure 4. Capped relative error (signed) for the BLYP functional for both BJ and CHG damping before (-D3) and after (-D3M) the refitting for the z ≤ 1.0 subset of the validation data set plus the SSI, BBI, and Water2510 validation databases. Red circles represent overbound DFT-D3 energies, while blue circles represent underbound DFT-D3 energies. Values within ±20% are not shown for clarity.



ASSOCIATED CONTENT

S Supporting Information *

The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.jpclett.6b00780. Refitted -D3M and -D3M(BJ) parameters, benchmarking details, DFT-D performance per functional, DFT-D performance per database (PDF) Database geometries, benchmark energy, DFT energies, DFT-D3 energies (ZIP)

diagram for the refitted results suggest that even if more complex damping forms were utilized that take into account the types of interaction, little improvement could be made. Thus, increasingly complex damping forms may not be an avenue of research worth pursuing. As the refitting performed here is clearly geared toward intermolecular interactions, the CONF data sets were used to ensure that refitting does not significantly impact thermochemistry. Across all functionals and damping forms, the refitted parameters only worsen the CONF accuracy by 3% on the average. For the worst case scenario, B97-D3(BJ) compared to B97-D3M(BJ), the CONF MCURE worsens by 7% (29 to 36%, respectively). This demonstrates that, at least for the CONF data set, damping forms and values have a much smaller impact on intramolecular dispersion. Until this point, all DFT values were computed with the counterpoise correction (CP) in the QZVP basis set to eliminate basis set superposition error. In comparison to the CP-corrected QZVP results, MCURE for the validation set deviates only mildly when the basis is reduced to CP aDZ (0− 2%) or non-CP QZVP (1−5%) but deviates considerably at non-CP aDZ (17−32%). Therefore, the basis for CP computations must be of at least aDZ quality (and for nonCP calculations of QZVP quality) in order to avoid serious degradation in the accuracy due to basis set errors. In summary, we have demonstrated that the established “minima cross section” benchmark data sets do not fully explore the van der Waals well at short-range and do not have enough diversity to adequately describe all the interactions occurring between side chains in proteins, or those on the water dimer potential energy surface. Several current benchmark sets have been extended to satisfy the former requirement, and novel data sets such as SSI were utilized to satisfy the latter. Thus, an important byproduct of our investigation is a massive set (8299 structures) of CCSD(T)-level interaction energies that is balanced with respect to both the



AUTHOR INFORMATION

Corresponding Authors

*E-mail: *E-mail: *E-mail: *E-mail:

[email protected]. [email protected]. [email protected]. [email protected].

Notes

The authors declare no competing financial interest.



ACKNOWLEDGMENTS D.G.A.S and K.P are supported by the Donors of the American Chemical Society Petroleum Research Fund, the NSF CAREER Award CHE-1351978, and the Intramural Grants Program at Auburn University. L.A.B. and C.D.S. are supported by a grant provided by the United States National Science Foundation (Grant No. CHE-1300497). The Center for Computational Molecular Science and Technology is funded through a NSF CRIF award (Grant No. CHE-0946869) and by the Georgia Institute of Technology.



REFERENCES

(1) Burns, L. A.; Vazquez-Mayagoitia, A.; Sumpter, B. G.; Sherrill, C. D. Density-Functional Approaches to Noncovalent Interactions: A Comparison of Dispersion Corrections (DFT-D), Exchange-Hole Dipole Moment (XDM) Theory, and Specialized Functionals. J. Chem. Phys. 2011, 134, 084107. (2) Gráfová, L.; Pitoňaḱ , M.; Ř ezác,̌ J.; Hobza, P. Comparative Study of Selected Wave Function and Density Functional Methods for 2201

DOI: 10.1021/acs.jpclett.6b00780 J. Phys. Chem. Lett. 2016, 7, 2197−2203

Letter

The Journal of Physical Chemistry Letters Noncovalent Interaction Energy Calculations Using the Extended S22 Data Set. J. Chem. Theory Comput. 2010, 6, 2365−2376. (3) Smith, D. G. A.; Patkowski, K. Toward an Accurate Description of Methane Physisorption on Carbon Nanotubes. J. Phys. Chem. C 2014, 118, 544−550. (4) Smith, D. G. A.; Patkowski, K. Benchmarking the CO2 Adsorption Energy on Carbon Nanotubes. J. Phys. Chem. C 2015, 119, 4934−4948. (5) Sherrill, C. D.; Takatani, T.; Hohenstein, E. G. An Assessment of Theoretical Methods for Nonbonded Interactions: Comparison to Complete Basis Set Limit Coupled-Cluster Potential Energy Curves for the Benzene Dimer, the Methane Dimer, Benzene-Methane, and Benzene-H2S. J. Phys. Chem. A 2009, 113, 10146−10159. (6) Zhao, Y.; Truhlar, D. G. Design of Density Functionals That Are Broadly Accurate for Thermochemistry, Thermochemical Kinetics, and Nonbonded Interactions. J. Phys. Chem. A 2005, 109, 5656−5667. (7) Zhao, Y.; Truhlar, D. G. Benchmark Databases for Nonbonded Interactions and Their Use To Test Density Functional Theory. J. Chem. Theory Comput. 2005, 1, 415−432. (8) Jurečka, P.; Šponer, J.; Č erný, J.; Hobza, P. Benchmark Database of Accurate (MP2 and CCSD(T) Complete Basis Set Limit) Interaction Energies of Small Model Complexes, DNA Base Pairs, and Amino Acid Pairs. Phys. Chem. Chem. Phys. 2006, 8, 1985−1993. (9) Zhao, Y.; Schultz, N. E.; Truhlar, D. G. Design of Density Functionals by Combining the Method of Constraint Satisfaction with Parametrization for Thermochemistry, Thermochemical Kinetics, and Noncovalent Interactions. J. Chem. Theory Comput. 2006, 2, 364−382. (10) Zhao, Y.; Truhlar, D. G. The M06 Suite of Density Functionals for Main Group Thermochemistry, Thermochemical Kinetics, Noncovalent Interactions, Excited States, and Transition Elements: Two New Functionals and Systematic Testing of Four M06-Class Functionals and 12 Other Functionals. Theor. Chem. Acc. 2008, 120, 215−241. (11) Kristyán, S.; Pulay, P. Can (Semi)local Density Functional Theory Account for The London Dispersion Forces? Chem. Phys. Lett. 1994, 229, 175−180. (12) Thanthiriwatte, K. S.; Hohenstein, E. G.; Burns, L. A.; Sherrill, C. D. Assessment of the Performance of DFT and DFT-D Methods for Describing Distance Dependence of Hydrogen-Bonded Interactions. J. Chem. Theory Comput. 2011, 7, 88−96. (13) Ř ezác,̌ 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. (14) Faver, J. C.; Benson, M. L.; He, X.; Roberts, B. P.; Wang, B.; Marshall, M. S.; Kennedy, M. R.; Sherrill, C. D.; Merz, K. M., Jr. Formal Estimation of Errors in Computed Absolute Interaction Energies of Protein-Ligand Complexes. J. Chem. Theory Comput. 2011, 7, 790−797. (15) Ř ezác,̌ J.; Riley, K. E.; Hobza, P. Extensions of the S66 Data Set: More Accurate Interaction Energies and Angular-Displaced Nonequilibrium Geometries. J. Chem. Theory Comput. 2011, 7, 3466−3470. (16) 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. (17) Becke, A. D.; Johnson, E. R. A Unified Density-Functional Treatment of Dynamical, Nondynamical, and Dispersion Correlations. J. Chem. Phys. 2007, 127, 124108. (18) Tkatchenko, A.; Scheffler, M. Accurate Molecular Van Der Waals Interactions from Ground-State Electron Density and FreeAtom Reference Data. Phys. Rev. Lett. 2009, 102, 073005. (19) Grimme, S.; Hansen, A.; Brandenburg, J. G.; Bannwarth, C. Dispersion-Corrected Mean-Field Electronic Structure Methods. Chem. Rev. 2016, 116, 5105−5154. (20) Dion, M.; Rydberg, H.; Schröder, E.; Langreth, D. C.; Lundqvist, B. I. Van der Waals Density Functional for General Geometries. Phys. Rev. Lett. 2004, 92, 246401. (21) Vydrov, O. A.; Van Voorhis, T. Nonlocal Van der Waals Density Functional: The Simpler the Better. J. Chem. Phys. 2010, 133, 244103.

(22) Grimme, S. Semiempirical GGA-Type Density Functional Constructed with a Long-Range Dispersion Correction. J. Comput. Chem. 2006, 27, 1787−1799. (23) Vydrov, O. A.; Scuseria, G. E. Assessment of a Long-Range Corrected Hybrid Functional. J. Chem. Phys. 2006, 125, 234109. (24) Lechner, C.; Sax, A. F. Adhesive Forces Between Aromatic Molecules and Graphene. J. Phys. Chem. C 2014, 118, 20970−20981. (25) Mardirossian, N.; Head-Gordon, M. Exploring the Limit of Accuracy for Density Functionals Based on the Generalized Gradient Approximation: Local, Global Hybrid, and Range-Separated Hybrid Functionals With and Without Dispersion Corrections. J. Chem. Phys. 2014, 140, 18A527. (26) Mardirossian, N.; Head-Gordon, M. Mapping the Genome of Meta-Generalized Gradient Approximation Density Functionals: The Search for B97M-V. J. Chem. Phys. 2015, 142, 074111. (27) Grimme, S.; Ehrlich, S.; Goerigk, L. Effect of the Damping Function in Dispersion Corrected Density Functional Theory. J. Comput. Chem. 2011, 32, 1456−1465. (28) Marshall, M. S.; Burns, L. A.; Sherrill, C. D. Basis Set CCSD(T) : Best Convergence of the Coupled-cluster Correction, δMP2 Practices for Benchmarking Non-covalent Interactions and the Attendant Revision of the S22, NBC10, HBC6, and HSG Databases. J. Chem. Phys. 2011, 135, 194102. (29) Ř ezác, J.; Riley, K. E.; Hobza, P. Benchmark Calculations of Noncovalent Interactions of Halogenated Molecules. J. Chem. Theory Comput. 2012, 8, 4285−4292. (30) Parker, T. M.; Sherrill, C. D. Assessment of Empirical Models Versus High-Acuracy Ab Initio Methods for Nucleobase Stacking: Evaluating the Importance of Charge Penetration. J. Chem. Theory Comput. 2015, 11, 4197−4202. (31) Smith, D. G. A.; Patkowski, K. Interactions between Methane and Polycyclic Aromatic Hydrocarbons: A High Accuracy Benchmark Study. J. Chem. Theory Comput. 2013, 9, 370−389. (32) Li, S.; Smith, D. G. A.; Patkowski, K. An Accurate Benchmark Description of the Interactions between Carbon Dioxide and Polyheterocyclic Aromatic Compounds Containing Nitrogen. Phys. Chem. Chem. Phys. 2015, 17, 16560−16574. (33) Bukowski, R.; Szalewicz, K.; Groenenboom, G. C.; van der Avoird, A. Predictions of the Properties of Water from First Principles. Science 2007, 315, 1249−1252. (34) Bukowski, R.; Szalewicz, K.; Groenenboom, G. C.; van der Avoird, A. Polarizable Interaction Potential for Water from Coupled Cluster Calculations. I. Analysis of Dimer Potential Energy Surface. J. Chem. Phys. 2008, 128, 094313. (35) Gruzman, D.; Karton, A.; Martin, J. M. L. Performance of Ab Initio and Density Functional Methods for Conformational Equilibria of CnH2n+2 Alkane Isomers (n = 4 − 8). J. Phys. Chem. A 2009, 113, 11974−11983. (36) Wilke, J. J.; Lind, M. C.; Schaefer, H. F.; Császár, A. G.; Allen, W. D. Conformers of Gaseous Cysteine. J. Chem. Theory Comput. 2009, 5, 1511−1523. (37) Reha, D.; Valdes, H.; Vondrásě k, J.; Hobza, P.; Abu-Riziq, A.; Crews, B.; de Vries, M. S. Structure and IR Spectrum of PhenylalanylGlycyl-Glycine Tripetide in the Gas-Phase: IR/UV Experiments, Ab Initio Quantum Chemical Calculations, and Molecular Dynamic Simulations. Chem. - Eur. J. 2005, 11, 6803−6817. (38) Csonka, G. I.; French, A. D.; Johnson, G. P.; Stortz, C. A. Evaluation of Density Functionals and Basis Sets for Carbohydrates. J. Chem. Theory Comput. 2009, 5, 679−692. (39) Goerigk, L.; Grimme, S. A General Database for Main Group Thermochemistry, Kinetics, and Noncovalent Interactions − Assessment of Common and Reparameterized (meta-)GGA Density Functionals. J. Chem. Theory Comput. 2010, 6, 107−126. (40) Jeziorski, B.; Moszyński, R.; Szalewicz, K. Perturbation Theory Approach to Intermolecular Potential Energy Surfaces of van der Waals Complexes. Chem. Rev. 1994, 94, 1887−1930. (41) Hohenstein, E. G.; Sherrill, C. D. Wavefunction Methods for Noncovalent Interactions. WIREs Comput. Mol. Sci. 2012, 2, 304−326. 2202

DOI: 10.1021/acs.jpclett.6b00780 J. Phys. Chem. Lett. 2016, 7, 2197−2203

Letter

The Journal of Physical Chemistry Letters

and Density Fitting Approximations. J. Chem. Phys. 2003, 118, 8149− 8160. (65) Polly, R.; Werner, H.-J.; Manby, F. R.; Knowles, P. J. Fast Hartree-Fock Theory Using Local Density Fitting Approximations. Mol. Phys. 2004, 102, 2311−2321. (66) Parker, T. M.; Burns, L. A.; Parrish, R. M.; Ryno, A. G.; Sherrill, C. D. Levels of Symmetry Adapted Perturbation Theory (SAPT). I. Efficiency and Performance for Interaction Energies. J. Chem. Phys. 2014, 140, 094106. (67) Papajak, E.; Zheng, J.; Xu, X.; Leverentz, H. R.; Truhlar, D. G. Perspectives On Basis Sets Beautiful: Seasonal Plantings of Diffuse Basis Functions. J. Chem. Theory Comput. 2011, 7, 3027−3034.

(42) Singh, N. J.; Min, S. K.; Kim, D. Y.; Kim, K. S. Comprehensive Energy Analysis for Various Types of π-Interaction. J. Chem. Theory Comput. 2009, 5, 515−529. (43) Johnson, E. R.; Becke, A. D. A Post-Hartree-Fock Model of Intermolecular Interactions. J. Chem. Phys. 2005, 123, 024101. (44) Chai, J.; Head-Gordon, M. Long-Range Corrected Hybrid Density Functionals with Damped Atom-Atom Dispersion Corrections. Phys. Chem. Chem. Phys. 2008, 10, 6615−6620. (45) Henderson, T. M.; Janesko, B. G.; Scuseria, G. E. Generalized Gradient Approximation Model Exchange Holes for Range-Separated Hybrids. J. Chem. Phys. 2008, 128, 194105. (46) Steinmann, S. N.; Wodrich, M. D.; Corminboeuf, C. Overcoming Systematic DFT Errors for Hydrocarbon Reaction Energies. Theor. Chem. Acc. 2010, 127, 429−442. (47) Seth, M.; Ziegler, T.; Steinmetz, M.; Grimme, S. Modeling Transition Metal Reactions with Range-Separated Functionals. J. Chem. Theory Comput. 2013, 9, 2286−2299. (48) Werner, H.-J.; Knowles, P. J.; Knizia, G.; Manby, F. R.; Schütz, M. Molpro: A General-Purpose Quantum Chemistry Program Package. WIREs Comput. Mol. Sci. 2012, 2, 242−253. (49) Turney, J. M.; Simmonett, A. C.; Parrish, R. M.; Hohenstein, E. G.; Evangelista, F. A.; Fermann, J. T.; Mintz, B. J.; Burns, L. A.; Wilke, J. J.; Abrams, M. L.; et al. PSI4: An Open-Source Ab Initio Electronic Structure Program. WIREs Comput. Mol. Sci. 2012, 2, 556−565. (50) Grimme, S. Semiempirical Hybrid Density Functional with Perturbative Second-Order Correlation. J. Chem. Phys. 2006, 124, 034108. (51) Becke, A. D. Density-Functional Thermochemistry. 3. The Role of Exact Exchange. J. Chem. Phys. 1993, 98, 5648−5652. (52) Stephens, P. J.; Devlin, F. J.; Chabalowski, C. F.; Frisch, M. J. Ab-Initio Calculation of Vibrational Absorption and CircularDichroism Spectra Using Density-Functional Force-Fields. J. Phys. Chem. 1994, 98, 11623−11627. (53) Becke, A. D. Density-Functional Exchange-Energy Approximation with Correct Asymptotic Behavior. Phys. Rev. A: At., Mol., Opt. Phys. 1988, 38, 3098−3100. (54) Lee, C.; Yang, W.; Parr, R. G. Development of the Colle-Salvetti Correlation-Energy Formula into a Functional of the Electron Density. Phys. Rev. B: Condens. Matter Mater. Phys. 1988, 37, 785−789. (55) Perdew, J. P. Density-Functional Approximation for the Correlation Energy of the Inhomogeneous Electron Gas. Phys. Rev. B: Condens. Matter Mater. Phys. 1986, 33, 8822−8824. (56) Adamo, C.; Barone, V. Toward Reliable Density Functional Methods without Adjustable Parameters: The PBE0 Model. J. Chem. Phys. 1999, 110, 6158−6170. (57) Ernzerhof, M.; Scuseria, G. E. Assessment of the Perdew-BurkeErnzerhof Exchange-Correlation Functional. J. Chem. Phys. 1999, 110, 5029−5036. (58) Perdew, J. P.; Burke, K.; Ernzerhof, M. Generalized Gradient Approximation Made Simple. Phys. Rev. Lett. 1996, 77, 3865−3868. (59) Weigend, F.; Ahlrichs, R. Balanced Basis Sets of Split Valence, Triple Zeta Valence and Quadruple Zeta Valence Quality for H to Rn: Design and Assessment of Accuracy. Phys. Chem. Chem. Phys. 2005, 7, 3297−3305. (60) Dunning, T. H., Jr. Gaussian-Basis Sets for Use in Correlated Molecular Calculations. 1. The Atoms Boron through Neon and Hydrogen. J. Chem. Phys. 1989, 90, 1007−1023. (61) Kendall, R. A.; Dunning, T. H., Jr.; Harrison, R. J. Electron Affinities of the 1st-Row Atoms Revisited - Systematic Basis Sets and Wave Functions. J. Chem. Phys. 1992, 96, 6796−6806. (62) Boys, S. F.; Bernardi, F. Calculation of Small Molecular Interactions by Differences of Separate Total Energies - Some Procedures with Reduced Errors. Mol. Phys. 1970, 19, 553−566. (63) van Duijneveldt, F. B.; van Duijneveldt-van de Rijdt, J. G. C. M.; van Lenthe, J. H. State-of-the-Art in Counterpoise Theory. Chem. Rev. 1994, 94, 1873−1885. (64) Werner, H.; Manby, F. R.; Knowles, P. J. Fast Linear Scaling Second-Order Møller-Plesset Perturbation Theory (MP2) Using Local 2203

DOI: 10.1021/acs.jpclett.6b00780 J. Phys. Chem. Lett. 2016, 7, 2197−2203