My PhD thesis involved determining kinetic isotope effects (KIE) for aromatic electrophilic substitution reactions in an effort to learn more about the nature of the transition states involved.[cite]10.1039/p29750001209[/cite] I learnt relatively little, mostly because a transition state geometry is defined by 3N-6 variables (N = number of atoms) and its force constants by even more and you get only one or two measured KIE per reaction; a rather under-defined problem in terms of data! So I decided to spend a PostDoc learning how to invert the problem by computing the anticipated isotope effects using quantum mechanics and then comparing the predictions with measured KIE.[cite]10.1021/ja00486a013[/cite] Although such computation allows access to ALL possible isotope effects, the problem is still under-defined because of the lack of measured KIE to compare the predictions with. In 1995 Dan Singleton and Allen Thomas reported an elegant strategy to this very problem by proposing a remarkably simple method for obtaining KIE using natural isotopic abundances.[cite]10.1021/ja00141a030[/cite] It allows isotope effects to be measured for **all** the positions in one of the reactant molecules by running the reaction close to completion and then recovering unreacted reactant and measuring the changes in its isotope abundances using NMR. The method has since been widely applied[cite]10.1021/ja109686[/cite],[cite]10.1021/ja205674x[/cite] and improved.[cite]10.1038/nchembio.352[/cite] Here I explore how measured and calculated KIE can be reconciled.

The original example uses the Diels-Alder cycloaddition as an example, with the 2-methylbutadiene component being subjected to the isotopic abundance. Although comparison with calculation on related systems was done at the time[cite]10.1021/ja00100a037[/cite] the computational methods in use then did not allow effects such as solvation to be included. I thought it might be worth re-investigating this specific reaction using more modern methodology (ωB97XD/Def2-TZVPP/SCRF=xylenes), giving an opportunity for testing one key assumption made by Singleton and Allen, *viz* the use of an internal isotope reference for a site where the KIE is assumed to be exactly 1.000 (the 2-methyl group in this instance). This assumption made me recollect my post on how methyl groups might not be entirely passive by rotating (methyl “flags”) in the Diels-Alder reaction between *cis*-butene and 1,4-dimethylbutadiene. I had concluded that post by remarking that *Rotating methyl groups should be looked at more often as harbingers of interesting effects*, which in this context may mean that such rotations may not be entirely isotope agnostic.

To start, I note that the *endo* (closed shell, *i.e.* non-biradical; the wavefunction is STABLE to open shell solutions) transition state obtained for this reaction[cite]10.14469/ch/191299[/cite],[cite]10.14469/ch/191301[/cite] has a computed dipole moment of 6.1D, just verging into the region where solvation starts to make an impact. Perhaps the most important conclusion drawn from Singleton and Allen’s original article[cite]10.1021/ja00486a013[/cite] was that the presence of an apparently innocuous 2-methyl substituent is sufficient to render the reaction asynchronous, which means that one C-C bond forms faster than the other. They drew this conclusion from observing that the inverse deuterium isotope effect was larger at C1 than C4, the difference being well outside of their estimated errors. The calculations indicate that the two bonds have predicted lengths of 2.197 (to C1) and 2.294Å (to C4) at the transition state. This means that an asynchronicity as small as Δ0.1Å can be picked up in measured isotope effects!

The calculated activation free energy is 19.2 kcal/mol (0.044M), which is entirely reasonable for a reaction occurring slowly at room temperature. The barrier for the *exo* isomer is 21.0 kcal/mol, 1.8 kcal/mol higher in free energy.^{†} The measured isotope effects are shown below with the predicted values in brackets. The colour code is green=within the estimated experimental error, red=outside the error.

The following observations can be made:

- The internal isotope reference assumed as 1.000 is reasonable for carbon, but the “rotating methyl groups” resulting from hyper conjugation between the C-H groups and the π system do have a small effect resulting in a predicted KIE of 0.996 rather than the assumed 1.000. This will impact upon all the other measured values to some extent.
- All the predicted
^{‡}^{13}C isotope effects agree with experiment within the error estimated for the latter. The calculation also has its errors, of which the most obvious is that harmonic frequencies are used rather than the more correct anharmonic values. - The
^{2}H isotope effects show more deviation. This might be a combination of the assumption that the internal Me reference has no isotope effect coupled with the use of harmonic frequencies for the calculation. - Although the
^{2}H values differ somewhat beyond the experimental error, the E/Z effects are well reproduced by calculation. The inverse isotope effect for the (Z) configuration is significantly larger in magnitude than for the (E) form, as was indeed noted by Singleton and Thomas. - So too is the asymmetry induced by the methyl group. The inverse isotope effects are greater for the more completely formed bond (to C1) than for the lagging bond (to C4). They are indeed a very sensitive measure of reaction synchronicity.

The pretty good agreement between calculation and experiment provides considerable reassurance that the calculated properties of transition states can indeed be subjected to *reality checks using experiment*. Indeed, it takes little more than a day to compute a complete set of KIEs, far less than it takes to measure them. One could easily argue that such a computation should accompany measured KIE whenever possible.

^{‡}This gives me an opportunity to extol the virtues of effective RDM (research data management). The two DOIs for the data include files containing the full coordinates and force constant matrices for both reactant and TS. Using these, one can obtain frequencies for any isotopic substitution you might wish to make in <1 second each, and hence isotope effects not computed here. One option might be to *e.g.* invert the reactant from the 2-methylbutadiene to the maleic anhydride and hence compute the isotope effects expected on this species (not reported in the original article) or to monitor instead the product.[cite]10.1021/ja9636348[/cite]

^{†}The KIE have only subtle small differences to the *endo* isomer; too small to assign the stereochemistry with certainty.

#### Acknowledgments

This post has been cross-posted in PDF format at Authorea.

Tags: Allen Thomas, calculated activation free energy, Chemistry, Dan Singleton, Deuterium, Diels–Alder reaction, Isotope, Isotopes, Kinetic isotope effect, Nuclear physics, Physical organic chemistry, shell solutions

You probably don’t want to leave out referencing our work with Ken Houk on predicting the isotope effects for this Diels-Alder. 10.1021/ja9615278

Your way of comparing predicted and experimental KIEs is a little off, due to the detail that the experimental KIEs made the assumption that the methyl group KIE was negligible. One you allow for this, ωB97XD/Def2-TZVPP/SCRF=xylenes doesn’t really perform as well as many other DFTs, including amusingly B3LYP/6-31G* gas phase.

Anyone who heard one of my talks between 1996 and 2005 would have seen slides comparing predicted and experimental isotope effects for this reaction. Over the years we have used that as a test case for how DFT methods performed with cycloadditions and the effect of doing the calculations in differing ways. Solvent models and bigger basis sets and dispersion corrections don’t make much difference. Many DFTs do fine, but notably M06-2X makes the TS too synchronous (I have seen enough in other cases to believe this is a systematic error). CVT/SCT does not make much difference here, except that it improves the two harder “inside” H/D KIEs a little, but it is painful to get numerical convergence.

There are of course lots and lots of papers making comparisons of calculated and experimental KIEs, thinking off the top of my head of papers from Saunders, Schramm, Paneth, Hess, Schaad, Houk, and Wolfsberg. We have 28 JACS papers doing this in one way or another.

How are you calculating your KIEs? You want to include a tunneling correction.

Re: KIEs. They are obtained from the free energies computed for the relevant isotope and have no tunnelling correction. Tunnelling corrections are conventionally applied to primary isotope effects for H; I also appreciate that they have been shown to be important for some 13C effects.

Thanks Dan for the background to the calculation of KIE. I did include a very early reference, dating back to 1978 when I believe the first calculations based on QM TS models for a complex reaction (the Diels-Alder) were reported; there were others from that period I did not cite (e.g. 10.1021/ja00493a008). As you suggest, there have been many since.

Yes, the computed KIE are what might be called “true” KIE, since they do not depend on the assumption of an internal site unaffected by them, and since this assumption turned out not to be quite 1.00, a small correction could be applied.

Regarding choice of functional, one could endlessly try any of the other 200+ available to seek small improvement. My choice of ωB97XD was based on experiences calibrating it for activation free energies against CCSD(T) calculations of the same system (doi: 10.1021/10.1021/ma300803b). I would think that any errors could just as well arise from other assumptions, such as harmonic frequencies, as from the choice of functional. Clearly, the B3LYP gas phase result may be one of those “right answer for the wrong reasons” cases that are quite common.

It would indeed also be good to see if anharmonic corrections affected calculated KIE significantly. I suspect not for 13C, but it might change the 2H ones. Unfortunately cubic force constants are far more arduous to compute for the models here than quadratic ones (indeed effectively impossible). It could be (perhaps has been) done for smaller reactions? Perhaps a suitable one for anharmonic study might be the HCl elimination reaction from chloroalkanes that I reported calculated KIE for back in 1981, doi: 10.1039/C39810000939. That has the advantage that the experiments are done in the gas phase.

A few points:

Calculating KIEs from the Gaussian free energies has a few disadvantages. The alternative is to calculate the KIEs from the Bigeleisen equation using programs that do so directly, such as quiver. It has been shown (I can find the ref if you want) that the former is more susceptible to errors in the low frequencies than the latter. Both benefit from cancellation of errors, so errors from the free-energy procedure are usually negligible, but the Bigeleisen equation is set up in a way that cancels errors in small numbers instead of errors in large numbers. A minor issue is that Gaussian does not quite give enough digits for 13C KIEs – a millionth of a hartree starts to make a difference. A third issue is that it is more work. From ordinary Gaussian outputs including a frequency calculation, I can have a complete set up isotope effects including tunneling corrections in about a minute.

To not include a tunneling correction for 13C KIEs is to ignore the literature that shows that tunneling corrections improve predictions. The effects are not negligible.

When a calculational prediction is wrong, it seems to be absolutely routine for the difference to be ascribed to some intractable or incalculable extra factor, rather than a direct error in calculated energies, geometries, or mechanism. This is done every day with entropy. Your worry about anharmonic corrections seems to fall into this category. Before you go down that road, it would seem that you should first explain a literature full of correct isotope effect predictions.

Finally, there is a rather important pattern that emerges when you do a lot of differing calculations on a reaction and check to see which correctly predict the KIEs. The calculations that give correct predictions also have about the same TS geometry. The calculations that give wrong predictions give differing geometries. The geometries might be “right for the wrong reason” but one can usually live with that. Errors are not due to the bogeyman, they are due to errors in transition states. All of this (with some caveats) is buried in that obscure JACS journal somewhere.

Thanks Dan. That is all useful information. To which information I would add the article by Wolfsberg, Saunders, Anet and Kronja in which the program THERMISTP (the successor to Quiver) is described; 10.1021/ja072375r. It appears the program can be requested by contacting its authors directly and this article also serves as an excellent tutorial in its use!