Archive for the 'DFT' Category

Dispersion in organic chemistry – a review and another example

The role of dispersion in organic chemistry has been slowly recognized as being quite critical in a variety of systems. I have blogged on this subject many times, discussing new methods for properly treating dispersion within quantum computations along with a variety of molecular systems where dispersion plays a critical role. Schreiner1 has recently published a very nice review of molecular systems where dispersion is a key component towards understanding structure and/or properties.

In a similar vein, Wegner and coworkers have examined the Z to E transition of azobenzene systems (1a-g2a-g) using both experiment and computation.2 They excited the azobenzenes to the Z conformation and then monitored the rate for conversion to the E conformation. In addition they optimized the geometries of the two conformers and the transition state for their interconversion at both B3LYP/6-311G(d,p) and B3LYP-D3/6-311G(d,p). The optimized structure of the t-butyl-substituted system is shown in Figure 1.

a: R=H; b: R=tBu; c: R=Me; d: R=iPr; e: R=Cyclohexyl; f: R=Adamantyl; g: R=Ph




Figure 1. B3LYP-D3/6-311G(d,p) optimized geometries of 1a, 2a, and the TS connecting them.

The experiment finds that the largest activation barriers are for the adamantly 1f and t-butyl 1b azobenzenes, while the lowest barriers are for the parent 1a and methylated 1c azobenzenes.

The trends in these barriers are not reproduced at B3LYP but are reproduced at B3LYP-D3. This suggests that dispersion is playing a role. In the Z conformations, the two phenyl groups are close together, and if appropriately substituted with bulky substituents, contrary to what might be traditionally thought, the steric bulk does not destabilize the Z form but actually serves to increase the dispersion stabilization between these groups. This leads to a higher barrier for conversion from the Z conformer to the E conformer with increasing steric bulk.


(1) Wagner, J. P.; Schreiner, P. R. "London Dispersion in Molecular Chemistry—Reconsidering Steric Effects," Angew. Chem. Int. Ed. 2015, 54, 12274-12296, DOI: 10.1002/anie.201503476.

(2) Schweighauser, L.; Strauss, M. A.; Bellotto, S.; Wegner, H. A. "Attraction or Repulsion? London Dispersion Forces Control Azobenzene Switches," Angew. Chem. Int. Ed. 2015, 54, 13436-13439, DOI: 10.1002/anie.201506126.


1b: InChI=1S/C28H42N2/c1-25(2,3)19-13-20(26(4,5)6)16-23(15-19)29-30-24-17-21(27(7,8)9)14-22(18-24)28(10,11)12/h13-18H,1-12H3/b30-29-

2b: InChI=1S/C28H42N2/c1-25(2,3)19-13-20(26(4,5)6)16-23(15-19)29-30-24-17-21(27(7,8)9)14-22(18-24)28(10,11)12/h13-18H,1-12H3/b30-29+

DFT &Schreiner Steven Bachrach 04 Jan 2016 No Comments

Benchmarking π-conjugation

With the proliferation of density functionals, selecting the functional to use in your particular application requires some care. That is why there have been quite a number of benchmark studies (see these posts for some examples). Yu and Karton have now added to our benchmark catalog with a study of π-conjugation.1

They looked at a set of 60 reactions which involve a reactant with π-conjugation and a product which lacks conjugation. A few examples, showing examples involving linear and cyclic systems, are shown in Scheme 1.

Scheme 1.

The reaction energies were evaluated at W2-F12, which should have an error of a fraction of a kcal mol-1. Three of the reactions can be compared with experimental values, and difference in the experimental and computed values are well within the error bars of the experiment. It is too bad that the authors did not also examine 1,3-cyclohexadiene → 1,4-cyclohexadiene, a reaction that is both of broader interest than many of the ones included in the test set and can also be compared with experiment.

These 60 reactions were then evaluated with a slew of functionals from every rung of Jacob’s ladder. The highlights of this benchmark study are that most GGA and meta-GGA and hybrid functionals (like B3LYP) have errors that exceed chemical accuracy (about 1 kcal mol-1). However, the range-separated functionals give very good energies, including ωB97X-D. The best results are provided with double hybrid functionals. Lastly, the D3 dispersion correction does generally improve energies by 10-20%. On the wavefunction side, SCS-MPs gives excellent results, and may be one of the best choices when considering computational resources.


(1) Yu, L.-J.; Karton, A. "Assessment of theoretical procedures for a diverse set of isomerization reactions involving double-bond migration in conjugated dienes," Chem. Phys. 2014, 441, 166-177, DOI: 10.1016/j.chemphys.2014.07.015.

DFT Steven Bachrach 14 Jan 2015 2 Comments

Becke’s Perspective on DFT

The Journal of Chemical Physics has produced a Special Topics issue on Advances in Density Functional Theory. I want to call to your attention the Perspective article by Becke titled “Perspective: Fifty years of density-functional theory in chemical physics”.1 Becke writes a personal account of the history of DFT and makes a number of interesting points and observations. He rightly notes that DFT is exact and we should more properly refer to our actual implementations as Density Functional Approximations (DFA). He also notes that use of the term ab initio as a synonym for wavefunction theory is inappropriate as DFT is just as ab initio as HF and post-HF theories.

A common perception about DFT (well, DFA) is that there is no way to systematically improve functionals. Becke exposes a true underlying logic that has driven much of DFA development.

Lastly, Becke is discouraged by the more recent developments that have included virtual orbitals, such as double hybrid methods. His approach is that true DFT is occupied orbitals only (for which he pointedly does not want to adopt the acronym OOO), and that developments that include the virtual orbitals might toll the “death knell” for DFT.

For those interested in a pretty accessible account of the history of DFT, Becke’s Perspective is an excellent place to get started.


(1) Becke, A. D. "Perspective: Fifty years of density-functional theory in chemical physics," J. Chem. Phys. 2014, 140, 18A301 DOI: 10.1063/1.4869598.

DFT Steven Bachrach 24 Sep 2014 No Comments

Computing OR: norbornenone

Optical activity is a major tool for identifying enantiomers. With recent developments in computational techniques, it is hoped that experiments combined with computations will be a powerful tool for determining absolute configuration. The recent work of Lahiri, et al. demonstrates the scope of theoretical work that is still needed to really make this approach broadly applicable.1

They prepared (1R,4R)-norbornenone 1 and measured its optical rotation in the gas phase and in dilute solutions of acetonitrile and cyclohexane. The specific rotations at three different wavelengths are listed in Table 1. Of first note is that though there is some small differences in solution, as expected, there really is substantial differences between the gas- and solution phases. Thus cautionary point 1: be very careful of comparing solution phase experimental optical activity with computed gas phase predictions.


Table 1. Experimental and computed specific rotation of 1.


355.0 nm

589.3 nm

633.0 nm

Gas phase













Acetonitrile solution









Cyclohexane solution









For the computations, the geometry of 1 was optimized at B3LYP/aug-cc-pVTZ (see Figure 1. The OR was computed at B3LYP with different basis sets, finding that the difference in the predicted specific rotation at 598.3nm differs only quite little (90.6 deg dm-1 (g/mL)-1) between the computations using aug-cc-pVTZ and aug-cc-pVQZ) suggesting that the basis set limit has been reached. The gas –phase computed values at B3LYP and CCSD are also shown in Table 1. Though these computations do get the correct sign of the rotation, they are appreciably off of the experimental values, and the percent error varies with wavelength. Cautionary point 2: it is not obvious what is the appropriate computational method to compute OR, and beware of values that seem reasonable at one wavelength – this does not predict a good agreement at a different wavelength.

Figure 1. Optimized geometry of 1 at B3LYP/aug-cc-pVTZ.

Lastly, computed solution values of the OR were performed with PCM and B3LYP. These values are given in Table 1. Again the agreement is poor. So cautionary point 3: computed (PCM) solution OR
may be in quite poor agreement with experiment.

Often the culprit to poor agreement between computed and experimental OR is attributed to omitted vibrational effects, but in this case, because 1 is so rigid, one might not expect too much error to be caused by the effects of vibrations. So the overall result – we need considerable work towards addressing how to accurately compute optical activity!


(1) Lahiri, P.; Wiberg, K. B.; Vaccaro, P. H.; Caricato, M.; Crawford, T. D. "Large Solvation Effect in the Optical Rotatory Dispersion of Norbornenone," Angew. Chem. Int. Ed. 2014, 53, 1386-1389, DOI: 10.1002/anie.201306339.


1: InChI=1S/C7H8O/c8-7-4-5-1-2-6(7)3-5/h1-2,5-6H,3-4H2/t5-,6+/m1/s1

DFT &Optical Rotation Steven Bachrach 25 Feb 2014 2 Comments

ORD of methyloxirane

Computing the optical rotation of simple organic molecules can be a real challenge. One of the classic problems is methyloxirane. DFT typically gets the wrong sign, let alone the wrong value. Cappelli and Barone1 have developed a QM/MM procedure where methyloxirane is treated with DFT (B3LYP/aug-cc-pVDZ or CAM-B3LYP/aubg-cc-pVDZ). Then 2000 arrangements of water about methyloxirane were obtained from an MD simulation. For each of these configurations, a supermolecule containing methyloxirane and all water molecules with 16 Å was identified. The waters of the supermolecule were treated as a polarized force field. This supermolecule is embedded into bulk water employing a conductor-polarizable continuum model (C-PCM). Lastly, inclusion of vibrational effects, and averaging over the 2000 configurations, gives a predicted optical rotation at 589 nm that is of the correct sign (which is not accomplished with a gas phase or simple PCM computation) and is within 10% of the correct value. The full experimental ORD spectrum is also quite nicely matched using this theoretical approach.


(1) Lipparini, F.; Egidi, F.; Cappelli, C.; Barone, V. "The Optical Rotation of Methyloxirane in Aqueous Solution: A Never Ending Story?," J. Chem. Theor. Comput. 2013, 9, 1880-1884, DOI: 10.1021/ct400061z.



DFT &Optical Rotation Steven Bachrach 15 May 2013 2 Comments

Benchmarking conformations: melatonin

Conformational analysis is one of the tasks that computation chemistry is typically quite adept at and computational chemistry is frequently employed for this purpose. Thus, benchmarking methods for their ability to predict accurate conformation energies is quite important. Martin has done this for alkanes1 (see this post), and now he has looked at a molecule that contains weak intramolecular hydrogen bonds. He examined 52 conformations of melatonin 1.2 The structures of the two lowest energy conformations are shown in Figure 1.




Figure 1. Structures of the two lowest energy conformers of 1 at SCS-MP2/cc-pVTZ.

The benchmark (i.e. accurate) relative energies of these conformers were obtained at MP2-F12/cc-pVTZ-F12 with a correction for the role of triples: (ECCSD(T)/cc-pVTZ)-E(MP2/cc-pVTZ)). The energies of the conformers were computed with a broad variety of basis sets and quantum methodologies. The root mean square deviation from the benchmark energies is used as a measure of the utility of these alternate methodologies. Of particular note is that HF predicts the wrong ordering of the two lowest energy isomers, as do some DFT methods that use small basis sets and do not incorporate dispersion.

In fact, other than the M06 family or double hybrid functionals, all of the functionals examined here (PBE. BLYP, PBE0, B3LYP, TPSS0 and TPSS) have RMSD values greater than 1 kcal mol-1. However, inclusion of a dispersion correction, Grimme’s D2 or D3 variety or the Vydrov-van Voorhis (VV10) non-local correction (see this post for a review of dispersion corrections), reduces the error substantially. Among the best performing functionals are B2GP-PLYP-D3, TPSS0-D3, DSD-BLYP and M06-2x. They also find the MP2.5 method to be a practical ab initio alternative. One decidedly unfortunate result is that large basis sets are needed; DZ basis sets are simply unacceptable, and truly accurate performance requires a QZ basis set.


(1) 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, DOI: 10.1021/jp903640h.

(2) Fogueri, U. R.; Kozuch, S.; Karton, A.; Martin, J. M. L. "The Melatonin Conformer Space: Benchmark and Assessment of Wave Function and DFT Methods for a Paradigmatic Biological and
Pharmacological Molecule," J. Phys. Chem. A 2013, 117, 2269-2277, DOI: 10.1021/jp312644t.


1: InChI=1S/C13H16N2O2/c1-9(16)14-6-5-10-8-15-13-4-3-11(17-2)7-12(10)13/h3-4,7-8,15H,5-6H2,1-2H3,(H,14,16)

DFT &MP Steven Bachrach 11 Apr 2013 2 Comments

Large water clusters and DFT performance

Truhlar has made a comparison of binding energies and relative energies of five (H2O)16 clusters.1 While technically not organic chemistry, this paper is of interest to the readership of this blog as it compares a very large collection of density functionals on a problem that involves extensive hydrogen bonding, a problem of interest to computational organic chemists.

The CCSD(T)/aug-cc-pVTZ//MP2/aug-cc-pVTZ energies of clusters 1-5 (shown in Figure 1) were obtained by Yoo.2 These clusters are notable not just for their size but also that they involve multiple water molecules involved in four hydrogen bonds. Truhlar has used these geometries to compute the energies using 73 different density functionals with the jun-cc-pVTZ basis set (see this post for a definition of the ‘jun’ basis sets). Binding energies (relative to 16 isolated water molecules) were computed along with the 10 relative energies amongst the 5 different clusters. Combining the results of both types of energies, Truhlar finds that the best overall performance relative to CCSD(T) is obtained with ωB97X-D, a hybrid GGA method with a dispersion correction. The next two best performing functionals are LC-ωPBE-D3 and M05-2x. The best non-hybrid performance is with revPBE-D3 and B97-D.

1 (0.0)

2 (0.25)

3 (0.42)

4 (0.51)

5 (0.54)

Figure 1. MP2/aug-cc-pVTZ optimized geometries and relative CCSD(T) energies (kcal mol-1) of (water)16 clusters 1-5. (Don’t forget to click on any of these molecules above to launch Jmol to interactively view the 3-D structure. This feature is true for all molecular structures displayed in all of my blog posts.)

While this study can help guide selection of a functional, two words of caution. First, Truhlar notes that the best performing methods for the five (H2O)16 clusters do not do a particularly great job in getting the binding and relative energies of water hexamers, suggesting that no single functional really stands out as best. Second, a better study would also involve geometry optimization using that particular functional. Since this was not done, one can garner little here about what method might be best for use in a typical study where a geometry optimization must also be carried out.


(1) Leverentz, H. R.; Qi, H. W.; Truhlar, D. G. "Assessing the Accuracy of Density
Functional and Semiempirical Wave Function Methods for Water Nanoparticles: Comparing Binding and Relative Energies of (H2O)16 and (H2O)17 to CCSD(T) Results," J. Chem. Theor. Comput. 2013, ASAP, DOI: 10.1021/ct300848z.

(2) Yoo, S.; Aprà, E.; Zeng, X. C.; Xantheas, S. S. "High-Level Ab Initio Electronic Structure Calculations of Water Clusters (H2O)16 and (H2O)17: A New Global Minimum for (H2O)16," J. Phys. Chem. Lett. 2010, 1, 3122-3127, DOI: 10.1021/jz101245s.

DFT &Hydrogen bond &Truhlar Steven Bachrach 25 Feb 2013 1 Comment

Benchmarked Dispersion corrected DFT and SM12

This is a short post mainly to bring to the reader’s attention a couple of recent JCTC papers.

The first is a benchmark study by Hujo and Grimme of the geometries produced by DFT computations that are corrected for dispersion.1 They use the S22 and S66 test sets that span a range of compounds expressing weak interactions. Of particular note is that the B3LYP-D3 method provided the best geometries, suggesting that this much (and justly) maligned functional can be significantly improved with just the simple D3 fix.

The second paper entails the description of Truhlar and Cramer’s latest iteration on their solvation model, namely SM12.2 The main change here is the use of Hirshfeld-based charges, which comprise their Charge Model 5 (CM5). The training set used to obtain the needed parameters is much larger than with previous versions and allows for treating a very broad set of solvents. Performance of the model is excellent.


(1) Hujo, W.; Grimme, S. "Performance of Non-Local and Atom-Pairwise Dispersion Corrections to DFT for Structural Parameters of Molecules with Noncovalent Interactions," J. Chem. Theor. Comput. 2013, 9, 308-315, DOI: 10.1021/ct300813c

(2) Marenich, A. V.; Cramer, C. J.; Truhlar, D. G. "Generalized Born Solvation Model SM12," J. Chem. Theor. Comput. 2013, 9, 609-620, DOI: 10.1021/ct300900e

Cramer &DFT &Grimme &Solvation &Truhlar Steven Bachrach 14 Jan 2013 No Comments

Monosaccharides benchmark

A comprehensive evaluation of how different computational methods perform in predicting the energies of monosaccharides comes to some very interesting conclusions. Sameera and Pantazis1 have examined the eight different aldohexoses (allose, alltrose, glucose, mannose, gulose, idose, galactose and talose), specifically looking at different rotomers of the hydroxymethyl group, α- vs. β-anomers, pyranose vs. furanose isomers, ring conformations (1C4 vs skew boat forms), and ring vs. open chain isomers. In total, 58 different structures were examined. The benchmark computations are CCSD(T)/CBS single point energies using the SCS-MP2/def2-TZVPP optimized geometries. The RMS deviation from these benchmark energies for some of the many different methods examined are listed in Table 1.

Table 1. Average RMS errors (kJ mol-1) of the 58 different monosaccharide structures for
different computational methods.


average RMS error























Perhaps the most interesting take-home message is that CEPA, MP2, the double hybrid methods and M06-2x all do a very good job at evaluating the energies of the carbohydrates. Given the significant computational advantage of M06-2x over these other methods, this seems to be the functional of choice! The poorer performance of the DFT methods over the ab initio methods is primarily in the relative energies of the open-chain isomers, where errors can be on the order of 10-20 kJ mol-1 with most of the functionals; even the best overall methods (M06-2x and the double hybrids) have errors in the relative energies of the open-chain isomers of 7 kJ mol-1. This might be an area of further functional development to probe better treatment of the open-chain aldehydes vs. the ring hemiacetals.


(1) Sameera, W. M. C.; Pantazis, D. A. "A Hierarchy of Methods for the Energetically Accurate Modeling of Isomerism in Monosaccharides," J. Chem. Theory Comput. 2012, 8, 2630-2645, DOI:10.1021/ct3002305

DFT &sugars Steven Bachrach 28 Nov 2012 No Comments

DSD-DFT – a double hybrid variation

I just returned from the Southwest Theoretical Chemistry Conference held at Texas A&M University. My thanks again to Steven Wheeler for the invitation to speak at the meeting and for putting together a very fine program and conference.

Among the many interesting talks was one by Sebastian Kozuch who reported on an interesting double hybrid methodology.1,2 Working with Jan Martin, they defined a procedure that Kozuch referred to as “putting Stefan Grimme into a blender”. They extend the double hybrid concept first suggested by Grimme that adds on an MP2-like correction functional. Kozuch and Martin substitute a spin-component scaled MP2 (SCS-MP2) model for the original MP2 correction. SCS-MP2 was also proposed by Grimme. Lastly, they add on a dispersion correction, an idea championed by Grimme too. The exchange-correlation term is defined as

EXC = cXEX DFT + (1 – cx)ExHF + cCECDFT + cOEOMP2 + cSESMP2 + s6ED

where cX is the coefficient for the amount of DFT exchange, cC the amount of DFT correlation, cC and cS the amount of opposite- and same-spin MP2, and s6 the amount of dispersion. They name this procedure DSD-DFT for Dispersion corrected, Spin-component scaled Double hybrid DFT.

In their second paper on this subject, they propose the use of the PBEP86 functional for the DFT components.2 Benchmarking against a variety of standard databases, including kinetic data, thermodynamic data, along with inorganic and weakly interacting systems, this method delivers the lowest mean error among a small set of functionals. Kozuch reported at the conference on a number of other combinations and should have a publication soon suggesting an even better method. Importantly, these DSD-DFT computations can be run with most major quantum codes including Orca, Molpro, Q-Chem and Gaussian (with a series of IOP specifications).

While double hybrid methods don’t have quite the performance capabilities of regular DFT, density fitting procedures offer the possibility of a significant reduction in computational time. These DSD-DFT methods are certainly worthy of fuller explorations.


(1) Kozuch, S.; Gruzman, D.; Martin, J. M. L. "DSD-BLYP: A General Purpose Double Hybrid Density Functional Including Spin Component Scaling and Dispersion Correction," J. Phys. Chem. C, 2010, 114, 20801-20808,
DOI: 10.1021/jp1070852

(2) Kozuch, S.; Martin, J. M. L. "DSD-PBEP86: in search of the best double-hybrid DFT with spin-component scaled MP2 and dispersion corrections," Phys. Chem. Chem. Phys. 2011, 13, 20104-20107, DOI: 10.1039/C1CP22592H

DFT Steven Bachrach 30 Oct 2012 1 Comment

Next Page »