A case study of the mechanism of alcohol-mediated Morita Baylis-Hillman reactions. The importance of experimental observations.

The mechanism of the Morita Baylis-Hillman reaction has been heavily studied in the literature, and a long series of computational studies have defined complete theoretical energy profiles in these reactions. We employ here a combination of mechanistic probes, including the observation of intermediates, the independent generation and partitioning of intermediates, thermodynamic and kinetic measurements on the main reaction and side reactions, isotopic incorporation from solvent, and kinetic isotope effects, to define the mechanism and an experimental mechanistic free-energy profile for a prototypical Morita Baylis-Hillman reaction in methanol. The results are then used to critically evaluate the ability of computations to predict the mechanism. The most notable prediction of the many computational studies, that of a proton-shuttle pathway, is refuted in favor of a simple but computationally intractable acid-base mechanism. Computational predictions vary vastly, and it is not clear that any significant accurate information that was not already apparent from experiment could have been garnered from computations. With care, entropy calculations are only a minor contributor to the larger computational error, while literature entropy-correction processes lead to absurd free-energy predictions. The computations aid in interpreting observations but fail utterly as a replacement for experiment.

Introduction

For simple reactions involving only a single kinetic step, the “reaction mechanism” is in general completely defined by the structure of the transition state. This structure can be probed by the many kinetics-based tools of classical experimental chemistry, including the determination of rate laws, substituent effects, solvent effects, kinetic isotope effects (KIEs), and activation parameters. For a two-step reaction, mechanistic studies are intrinsically less decisive as the reaction now involves two transition states plus an intermediate. For mechanisms involving more steps, the complexities are multiplied. Often only one of the transition states, that for the rate-limiting step, can be scrutinized by kinetic probes, and intermediates are usually not directly observable. For many important multistep reactions, experimental studies can provide only limited glimpses of the mechanism.

For the understanding of complex reactions, the rise of computational mechanistic chemistry has arguably been the most important advance ever. The combination of reasonably accurate density functional theory (DFT) methods and ever-increasing computational power has stimulated the application of this technology on a broad front.[1] Few reactions, if any, are considered too complicated for computational study. Such studies then provide apparently complete mechanisms, including the geometries and energies of every intermediate and transition state. This level of detail is beyond the most ambitious dreams of classical experimental mechanistic chemistry.

This impressive accomplishment also constitutes a potential problem. That is, those mechanistic details that cannot be discerned from experimental studies are also not directly confirmable, or falsifiable, by experimental studies. The argument for the accuracy of such studies is usually an indirect one, most often based on the general accuracy of the potential energy surface for simpler problems or when compared with higher-level calculations. This scientific approach can go wrong on multiple levels. At one level, the accuracy of a theoretical method for some other problem may not imply accuracy for the problem at hand. At a second level, even a perfectly accurate potential energy surface may be quite misleading in comparison to the decisive free energy surface, and the allowance for entropy may be inaccurate or else impractical to achieve accurately. At a more human level, calculational studies do not speak to mechanistic possibilities that were not explored. In a complex reaction, possible mechanisms may easily be missed. Some conventional mechanistic steps, particularly proton transfers to or from solvent, are sufficiently intractable computationally that they may be ignored. Finally, the paradigms used to interpret computational mechanistic results, particularly statistical rate theories, may not be accurately applicable to a system under study, even for common organic reactions in solution.[2] It should be recognized that the goal of accuracy has been a central feature of computational mechanistic chemistry. No small effort has been exerted in this endeavor. However, it must also be recognized that the accuracy of many studies is ultimately both uncertain and unexamined.

We describe here a case study of a typical complex reaction, the alcohol-mediated Morita Baylis–Hillman (MBH) reaction,[3] using a full gamut of experimental mechanistic probes as well as a full standard computational study using the two most popular DFT methods augmented by high-level calculations on model reactions. The MBH mechanism in general outline (Scheme 1) is uncontroversial, consisting of the “addition” step by an activating nucleophile to afford 1, the carbon–carbon bond-forming “aldol” step affording 2, and an “elimination” step (either concerted or multistep) to afford the product.[4] However, the multicomponent nature of the MBH reaction and multistep nature of its mechanism provides fodder for many complications that affect experimental observations. We viewed the MBH mechanism as a special opportunity for mechanistic study because its individual steps are amenable to detailed scrutiny using many experimental probes, including the observation of intermediates, the independent generation and partitioning of intermediates, thermodynamic and kinetic measurements on the main reaction and side reactions, isotopic incorporation from solvent, and KIEs. The inferences from experimental studies are then compared with computational predictions. Our conclusions are pessimistic from one perspective; the computational studies are arguably more misleading than enlightening. It is not clear to us that any reliably accurate information that was not already apparent from experiment could have been garnered from calculations. The problem of calculating the entropy along the mechanistic pathway, however, does not appear as daunting as we originally expected. Our results highlight issues for careful consideration with regard to a broad genre of papers in the literature.

Scheme 1

The MBH reaction has been usefully catalyzed or promoted by tertiary amines, phosphines, oxygenated bases, Lewis acids, metals, water, alcohols, high pressure, ultrasound, autocatalysis, and even the use of lower temperatures.[5] This complexity underscores the importance of mechanistic understanding for the rational control of reactions and development of new reactions. There has thus been considerable interest in the MBH mechanism. Carbon–carbon bond-forming steps tend to have higher barriers than proton-transfer steps, so the aldol step might have been expected to be rate-limiting in the mechanism. This idea was supported by Hill and Isaacs for the DABCO-catalyzed reaction of acrylonitrile with acetaldehyde on the basis of third-order kinetics (rate = k[MeCHO][acrylonitrile][DABCO]), pressure dependence studies, and an H/D KIE for the α-position of acrylonitrile of 1.03 ± 0.1.[4] Third-order kinetics was also observed for reactions of acrylate esters with pyridinecarboxaldehydes.[6]

Much later, McQuade and co-workers stood the simplistic picture of the MBH mechanism on its head by the unexpected finding that the reaction of acrylates with aryl aldehydes in aprotic solvents was overall fourth order.[7] McQuade additionally observed large H/D KIEs for the α-position of acrylates. These observations were inconsistent with a rate-limiting aldol step. In its place, McQuade proposed that a rate-limiting elimination step was aided by a second molecule of aldehyde in a hemiacetal intermediate, as in 3.

The acceleration of the MBH reaction by alcohols and water has long been noted by many groups.[4b,6,8−10] Hill and Isaacs proposed that alcohols acted by hydrogen bonding that promoted the aldol step. Aggarwal and Lloyd-Jones observed that the reaction of methyl acrylate (MA) with benzaldehyde exhibited autocatalysis. From this, they proposed that the product alcohol was acting as a shuttle to transfer a proton from the α-position of 2 to the alkoxide via a six-membered cyclic transition state, as in 4.

A series of 11 papers from multiple groups has studied the MBH mechanism computationally.[11−19] Every paper that examined the issue, a total of seven, supported the Aggarwal/Lloyd-Jones proton shuttle depicted in 4, and this prediction was a highlight of most of these papers.[11,12,14a,15a,17] Large computational error is evident in some of these papers,[13b,14b,15b] but several of the groups undertook substantial and respected approaches to minimizing error. Sunoj chose his DFT method (MPW1K) on the basis of comparisons with high-level CBS-4M calculations in computational models.[17] Aggarwal and Harvey employed G3MP2 calculations on a model system to calibrate their B3LYP results.[11] Harvey later studied in detail the ability of diverse computational methods to predict the barrier for an MBH reaction, and he recognized explicitly the difficulty of predicting rate constants quantitatively.[18] Cantillo and Kappe chose M06-2X calculations on the basis of detailed experimental thermodynamics.[12]

The approaches to error minimization employed in these works are typical of the better computational mechanistic studies. They are clearly the result of recognition of the potential for computational error. However, the actual errors in the theoretical mechanism, free energies, enthalpies, and entropies along the reaction pathway were unknown for either this or any comparably complex organic mechanism. The multimolecular nature of the MBH mechanism made the ability of calculations to predict entropy changes a particular concern. Our studies were initiated with the goal of remedying the general ignorance of error for a specific example of an MBH reaction. In this way we sought to gain insight into error in the broader perspective of computational mechanistic studies.

Results and Discussion

The prototypical MBH reaction of p-nitrobenzaldehyde (5) with MA in methanol catalyzed by DABCO (6) was chosen for study. This reaction cleanly affords the adduct 7 in both methanol and DMSO at 25 °C, with the DMSO reaction requiring extended reaction times. The kinetics for the DMSO reaction as well as its α-position and aldehydic H/D KIEs had been studied by McQuade.[7]

MBH Thermodynamics

The reversibility of the MBH reaction enforces substantial limitations on the reaction scope. The reaction of 5 with MA is relatively exergonic among MBH reactions but it does not proceed to completion. Cantillo and Kappe had studied the equilibrium of 5/MA with 7 in methanol at a range of temperatures to obtain thermodynamic parameters,[12] but there were unrecognized complications in their study. One problem is that 5 in methanol is in rapid equilibrium with its hemiacetal (p-O2NC6H4CH(OH)(OMe), present at 79% at 25 °C). A second problem is that slow side reactions occur at elevated temperatures. By NMR analysis, we were able to measure directly the concentrations of 5, MA, and 7 in equilibrium mixtures obtained from both the forward reactions of 5 and MA and the reverse reaction of purified 7 in d4-methanol using 30 mol % of DABCO. The equilibrium constants at 22 and 60 °C were 860 M–1 and 66 M–1, respectively. Our data fit with ΔH° = −13.2 kcal/mol and ΔS° = −31 e.u., putting the ΔG° = −3.9 kcal/mol at 25 °C.

The Addition Step: Shunt Processes and Thermodynamics

The addition step is normally depicted simplistically as in Scheme 1, but a complication is that the zwitterionic 1, or more specifically 8, can be protonated by alcohols or water. The reversible protonic equilibrium with methanol forming 10 involves proton transfers between heteroatoms and should be rapid,[20] but this O-protonation is hidden from experimental detection. The more interesting and experimentally tractable process is the C-protonation affording 9.

To probe the formation of 9 under the reaction conditions, the reaction in d4-methanol was followed by 1H NMR. Deuterium incorporation into the unreacted MA was extensive; by the time that the formation of 7 was 18% complete, 85% of the MA was deuterated. Isotopic exchange in the absence of DABCO is negligible. This shows that the C-protonation of 8 is faster than product formation. At a series of points early in the reaction (Figure 1), the deuterium incorporation was an approximate factor of 5 greater than the formation of 7. If it is assumed that molecules of 9 with one D and one H in the α-position most often lose H in returning to 8, then the factor of 5 represents the relative rate for the C-protonation of 8 versus the rate-limiting step for product formation. Since the subsequent aldol and elimination steps should be normal, barriered processes, the C-protonation of 8 must also have a substantial barrier. This is as expected for a proton transfer to carbon, but the point is pertinent to later discussion.

Figure 1

Example kinetics runs, showing reactions of 5 with MA or α-d-MA in methanol or d4-methanol. The marked points are for experimental observations. The solid lines are theoretical curves based on the rate law rate = k[5][MA][DABCO], with k being the value listed, derived by fitting to the experimental points. The green solid line represents a fit to the initial four points for the reaction of MA in d4-methanol; later points fall off the curve due to deuterium incorporation in the MA α-position.

While 8 is too unstable to be observed, the cationic adduct 9 is readily observable in the reaction of the hydrochloride salt of DABCO with MA (catalyzed by DABCO free base). This allows us to assess the stability of 8 using the thermodynamic cycle of Scheme 2. In this cycle, the unobservable equilibrium of MA and DABCO with 8 is related to the observable equilibrium of MA and DABCO-H+ simply by the difference in the acidity of 9 versus DABCO-H+.

Scheme 2

Thermodynamic Cycle Defining the Stability of 8

Equilibrium constants KNH for the conversion of MA/DABCO-H+ to 9 in d4-methanol were determined at a series of temperatures by NMR. A complication in this observation was that the equilibration was too slow to carry out within the spectrometer. Instead, samples were rapidly warmed or cooled to ambient temperature and a series of spectra were taken to allow extrapolation of the concentrations back to the original mixture. At temperatures of 0, 22, 40, and 60 °C, the KNH values were 1170, 260, 94, and 34.3 M–1, respectively. A plot of ln KNH versus 1/T gave ΔH° = −10.6 ± 0.3 kcal/mol and ΔS° = −24.9 ± 0.9 e.u.

The calculation of the thermodynamics of formation of 8 now requires an estimate of the difference in pKa of DABCO-H+ versus 9. This difference was assessed from the kinetic acidity of 9 when deprotonated by DABCO. Based on deuterium exchange into 9 in d4-methanol, ignoring any secondary KIE or internal return, and allowing for the two nitrogens in DABCO, the rate constant for deprotonation of 9 by DABCO was 7 × 10–4 M–1 s–1. To translate this rate constant into an equilibrium constant using Marcus theory,[21] we used the Guthrie equation,[22] log k = 10 – b[1 – ((log K)/4b)]2, with the parameter b set as 8.3 on the basis of Bernasconi’s observation of an intrinsic rate constant of 101.7 for the very similar deprotonation of a cationic ketone, 2-acetyl-1-methylpyridinium ion, by amines in 50% DMSO-H2O.[23] This leads to a log K (for deprotonation of 9 by DABCO) of −8.6. The pKa of DABCO is 8.8 in water,[24] and it changes little with solvent (9.06 in DMSO[25]) as is normal for cationic acids. This leads to a pKa ≈ 17.4 for 9. This value seems reasonable when it is considered that the nearby cationic charge in 2-acetyl-1-methylpyridinium ion lowers the pKa of a ketone by approximately 8.3 pKa units[23] and that the pKa of 17.4 is 8.2 less than the ethyl acetate pKa of 25.6.[26] An upper limit on log K for deprotonation of 9 by DABCO can be set by taking the reverse reaction as being diffusion controlled with a rate constant of approximately 1010 M–1 s–1. This would place log K at no more than −13.2, so the pKa of 9 can be no more than about 22. Since the protonation of the zwitterionic 8 by DABCO-H+ involves considerable reorganization, its rate is likely far less than diffusion controlled, bringing the pKa of 9 toward the first estimate. In support of this, it was noted above that the C-protonation of 8 under MBH conditions must have a substantial barrier to account for the similar rates of deuterium incorporation and MBH product formation. Another approach to assessing the acidity of 9 was a comparison of its rate of deuterium exchange with that of 3-pentanone under identical conditions. The rate constant for exchange into 3-pentanone was 180 times slower. The intrinsic barrier for deprotonation of simple ketones is lower than that for cationic ketones, so the difference in kinetic acidities suggests that 9 is several pKa units more acidic than 3-pentanone. The pKa of 3-pentanone in water is 19.9; it would be somewhat higher in methanol, but the pKa of the cationic 9 should be relatively solvent-independent. Overall, we will take log K for deprotonation of 9 by DABCO at −8.6, but allow that there is an uncertainty in this number of easily ±1. From this log K and the observed ΔG° for formation of 9, ΔG° for formation of 8 is approximately +8.6 kcal/mol.

Based on the rate constant of 8.2 × 10–4 M–1 s–1 at 25 °C for deuteration of MA under the MBH reaction conditions in d4-methanol (kobs = 1.4 × 10–4 s–1, [DABCO] = 0.17 M, rate = k[MA][DABCO]), the barrier for the formation of 9 is 21.7 kcal/mol. Unfortunately this barrier does not reflect the barrier for formation of 8, as the rate-limiting step in the formation of 9 is protonation of 8. This is known because DABCO-H+ catalyzes the formation of 9; the rate of formation of 9 in the presence of 0.667 M DABCO-H+ is approximately 100-fold faster than the rate of deuterium incorporation into MA, which requires the intermediacy of 9, in the absence of DABCO-H+.

A summary of the free-energy profile for the addition step of the MBH mechanism as derived from these observations is shown in Figure 2. In this profile, the standard state for methanol is taken as neat methanol while the standard state for all other compounds is 1 M. The free-energy barrier for formation of 8 remains unknown, though the thermodynamics for formation of 8 and the barrier for formation of 9 provide lower and upper bounds for this barrier. Under the MBH reaction conditions, the concentration of DABCO-H+ is approximately 10–4 M based on the autoprotolysis constant of methanol (10–16.71)[27] and an assumed pKa of DABCO-H+ in methanol of 9. Because the DABCO-H+-catalyzed process would be slow compared to the observed rate of deuterium incorporation, the formation of 9 must be mediated by the methanol and afford MeO–. The free energy of this combination is calculable from the methanol autoprotolysis constant and the pKa of DABCO-H+.

Figure 2

Experimental and calculated free energies for the addition step of the MBH mechanism. Energies are in kcal/mol and entropies are in e.u. The barriers for protonation of 8 were not calculated.

MBH Kinetics, Activation Parameters, and H/D KIEs

The progress of a series of MBH reactions of 5 with MA in methanol was followed by analysis of worked-up aliquots by NMR or HPLC. The observations were then kinetically simulated, and they fit well with a kinetic model that was first-order in 5 and first-order in MA (see Figure 1 and the Supporting Information (SI)). All relative rates were obtained from side-by-side reactions conducted identically.

A reaction with half the normal concentration of DABCO was slower by a factor of 2.0 ± 0.1, while a reaction with double the normal concentration of DABCO went faster by a factor of 1.9 ± 0.1. This indicates that the reaction is first order in DABCO. Addition of 30 mol % of the hydrochloride salt of DABCO had the effect of slowing the reaction by 30 ± 8% while adding 60 mol % slowed the reaction by 59 ± 8%. This result is as would be expected if most of the DABCO-H+ were rapidly converted to the less reactive 9 but otherwise the addition of the buffering DABCO-H+ had no effect of the rate. The rate was also unchanged in the presence of 30 mol % of Proton-sponge. These observations indicate that the number of protons in the rate-limiting transition state is the same as that in the starting materials. In other words, the total charge in the rate-limiting transition state is zero. All of our observations fit with rate = k[5][MA][DABCO]. The kinetics of course do not discern whether additional solvent molecules are specifically involved. The choice of neat methanol as its standard state avoids the need to adjust activation parameters for mechanisms involving additional methanol molecules.

A series of kinetics measurements were conducted with careful temperature control at temperatures ranging from −21.3 to 63.7 °C. A striking feature of these results is that the rate constant reaches a maximum near room temperature. An Eyring plot of the results (Figure 3) is decidedly and reproducibly nonlinear. As will be supported by later observations, the curvature in the Eyring plot is consistent with a reaction involving competitive rate-limiting steps having significantly different ΔS⧧. In such a case, the step with the more favorable (less negative) ΔS⧧ and higher ΔH⧧ would dominate the barrier at low temperatures while the step with the less favorable (more negative) ΔS⧧ and lower (negative in this case) ΔH⧧ would dominate the reaction at high temperatures, in accord with what is observed. If the reaction involved two separate, independent mechanisms with differing ΔS⧧, the curvature of the Eyring plot would be in the opposite direction!

Figure 3

Example Eyring plot based on kinetics runs from −21.3 to 63.7 °C. The solid line is simulated on the basis of eq 1

When two sequential steps are competitively rate-limiting and the steady-state approximation applies, the observed rate constant is governed by eq 1, where ΔG1⧧ and ΔG2⧧ are the total barriers versus starting materials for the two steps. The observed Eyring plot could then be simulated well with ΔH1⧧ = 12 ± 2 kcal/mol, ΔS1⧧ = −27 ± 9 e.u., ΔH2⧧ = −2.3 ± 1.3 kcal/mol, and ΔS2⧧ = −79 ± 5 e.u. Independent data gave similar results (see the SI). The errors in the fit ΔH⧧ and ΔS⧧ are notably not independent. For example, a positive error in ΔS1⧧ would imply a negative error in ΔS2⧧. The ΔG⧧ values are more precise. The ΔG1⧧ and ΔG2⧧ at 25 °C would be 20.2 ± 0.3 and 21.2 ± 0.2 kcal/mol, respectively.

When the reaction is conducted in d4-methanol, the kinetics depend on how the reaction is conducted. If the reaction is initiated by adding the DABCO last, the earliest part of the reaction involves unlabeled MA. The initial rate is then nearly equal to that of the reaction in unlabeled methanol (Figure 1), with the observed solvent kH/kD = 0.96 ± 0.1. The near-unity solvent KIE indicates that there is no proton transfer of hydroxylic protons in the majorly rate-limiting step. Qualitatively, this would appear to weigh strongly against a transition state of the type proposed by Aggarwal and Lloyd-Jones, as depicted in 4. It would also confute the calculational support for such a structure seen in the many computational studies. This interpretation will be considered in more detail below in the light of direct studies of the elimination reaction.

As the reaction in d4-methanol proceeds, it slows down as deuterium is incorporated into the α-position of MA. To determine the H/D KIE for the α-position of MA, the MA was first equilibrated with excess d4-methanol using DABCO in the absence of 5, then the reaction was initiated by adding 5. The resulting kH/kD was 3.1 ± 0.3, indicative of removal of an H/D from the α-position in the majorly rate-limiting step. From this last observation, we adopted the working hypothesis that ΔG2⧧ corresponds to the proton transfer of the elimination step while ΔG1⧧ corresponds to the aldol step. If true, it would be predicted that the rate-limiting step would change to the aldol step at low temperatures, where ΔG1⧧ dominates the barrier, and the H/D KIE should drop. This was found to be the case; at −20 °C the H/D KIE was only 1.1 ± 0.2. It will be seen that diverse other evidence supports the hypothesis.

13C KIEs

The 13C KIEs for the reaction of 5 with MA were determined at natural abundance by NMR methodology.[28] Duplicate independent reactions in both DMSO and methanol were taken to 77–80% conversion of 5. Reisolated aldehyde samples were then analyzed by 13C NMR in comparison with samples of the original aldehyde. The carbons meta to the aldehyde on the aromatic ring were treated as an internal standard with the assumption that their isotopic composition did not change. From the reaction conversions and the changes in the isotopic composition, the 13C KIEs were calculated as previously described.[28] Due to a long relaxation time and the width of its 13C peak, the para position in 5 could not practically be quantitated reliably.

The KIEs for 5 in DMSO and in methanol or d4-methanol are summarized in Figure 4a,b. In each case a significant but modest 13C KIE was observed for the aldehydic carbon. The remaining 13C KIEs observed were near unity, as would be expected for centers unchanged by the reaction. At 1.009, the carbonyl carbon KIE in DMSO is smaller than normally associated with a primary 13C KIE. The qualitative interpretation of this KIE is that the carbonyl carbon has undergone some process that has modified this center, but that no bond is being made or broken at this center at the transition state for the rate-limiting step. This is as would be expected for the McQuade mechanism. Because of concern over the effect of water or other hydroxylic impurities (including the product) on the reaction, two additional DMSO experiments were conducted, one taking precautions to minimize the presence of water, and a second with 1% water added to the reaction. The former had no impact on the KIEs while the latter led to a very slight increase in the KIE at the aldehydic carbon.

Figure 4

13C KIEs for the DABCO-catalyzed MBH reaction of 5 with MA at 25 °C. KIEs marked with * are for a reaction with precautions taken to minimize water, while KIEs marked with † are for a reaction with 1% water in DMSO used as solvent. KIEs marked ‡ and # were measured in d4-methanol, and the # signifies that the MA was pre-equilibrated with the d4-methanol before adding DABCO.

The KIEs in methanol and d4-methanol were more interesting. In methanol, the aldehydic carbon 13C KIE of 1.015–1.016 was significantly larger than it was in DMSO. However, it is still smaller than the large primary 13C KIE that would be expected if the addition to the aldehyde became rate-limiting. (See below for quantitative predictions.) If the DMSO and methanol reactions both had purely elimination processes as their rate-limiting step, then there would be no obvious explanation as to why the KIEs differ.

This line of reasoning supports the involvement of competitive rate-limiting steps in the methanol reaction. When one of the possible competitive rate-limiting steps involves proton transfer, a rigorous test for such kinetic complexity involves carrying out the reaction with a deuteron transferred instead of a proton, looking for a change in a carbon KIE.[29] The idea behind such a study is that the primary H/D KIE leads to a change in the relative importance of the mixed rate-limiting steps. In the event, this test was performed initially by simply carrying out the reaction in d4-methanol, and the aldehydic carbon 13C KIE decreased to 1.012. A flaw in this experiment was that a portion of the reaction occurred before high incorporation of deuterium into the MA. When the reaction was carried out with a pre-equilibration of the MA with d4-methanol in the presence of DABCO, the 13C KIE decreased to 1.009–1.011, which is indistinguishable from the DMSO KIEs. These observations strongly support the involvement of competitive rate-limiting steps.

The MA 13C KIEs in DMSO were straightforwardly measured. Reactions taken to ∼80% conversion were quenched by the addition of benzoic acid. The recovered unreacted MA was then analyzed by NMR in comparison with the original MA using the methyl carbon as the internal standard for quantitation. This process did not work for reactions in methanol due to problems with transesterification and recovery of the MA. As an alternative, samples of 7 from reactions taken to low conversion were analyzed versus samples taken to 100% conversion of the MA. Due to transesterification, the methyl carbon could not be used for quantitation but the negligible KIE in the aromatic carbons of 5 made the aromatic carbons of 7 suitable for use as internal standards. In d4-methanol, we unfortunately could not obtain MA KIEs due to the NMR complications associated with incorporation of deuterium.

The MA KIEs are summarized in Figure 4c,d. The β-carbon KIE is significantly inverse in each case. 13C KIEs of this magnitude suggest a pre-equilibrium converting the carbon to the more constraining potential energy well associated with sp3 hybridization, followed by a rate-limiting step that is unrelated to this carbon. If the elimination is rate-limiting, the inverse β-carbon KIE appears to exclude the concerted (E2) mechanism depicted in 3. Rather, the elimination would have to occur by a rate-limiting proton transfer followed by a faster loss of DABCO as a separate step in an E1cb(irr) process. The α-carbon KIE follows the pattern seen for the aldehydic carbon of 5: small in DMSO, and larger though still relatively small in methanol. These low KIEs are initially surprising since every reasonable mechanism involves some bonding change at the α-carbon in the rate-limiting step, but some insight into these KIEs will be obtained with the aid of calculations below.

The Aldol Step: Transition Structures, Predicted 13C KIEs, and Independent Experimental Energetics

For the aldol step, neither the starting 8 nor the product 11 or its protonated form 12 could be observed. This precludes direct experimental information about the step. However, a quantitative interpretation of the KIEs provides an independent if indirect estimate of the barrier for the aldol step. This required the aid of computations. Three computational approaches were explored. B3LYP[30] and M06-2X[31] calculations were carried out using a PCM solvent model[32] for methanol, while the M06-2X calculations were also performed using an SMD methanol solvent model.[33] Full optimization and a 6-31+G** basis set were used in all calculations unless otherwise stated.

A series of 12 aldol-step transition structures were located using each computational method. The 12 possibilities within each series arise from three rotational orientations of the aldehyde, attack on either the re or si face of the aldehyde, and the reaction of Z versus E isomers of 8. Some additional transition structures involving alternative orientations of the DABCO moiety were located but these were much higher in energy. The lowest-energy transition structures 13a–c (Figure 5) orient the DABCO moiety cis to the enolate oxygen and distal to the approaching aldehyde. Notably, the three computational approaches differ substantially in the geometry and early versus late character of the transition structures, and they disagree over the preferred diastereomer. As will be seen below, the B3LYP barrier is most accurate for this step but the M06-2X/PCM calculations provide better overall energetics, so the energetics provide no clear guidance as to which structure is better.

Figure 5

Lowest-energy transition structures and predicted 13C KIES at 25 °C for the reaction of 8 with 5.

The transition structures define the 13C KIEs that would be expected for a rate-limiting aldol step. The 13C KIEs predicted from conventional transition state theory were calculated from the scaled theoretical vibrational frequencies by the method of Bigeleisen and Mayer.[34] Tunneling corrections were applied using the one-dimensional infinite parabolic barrier model.[35] Such KIE predictions including a one-dimensional tunneling correction have proven highly accurate, so long as the calculation accurately depicts the mechanism and transition state geometry.[36]

The resulting predicted KIEs (Figure 5) are far from the experimental values. If the aldol step were fully rate-limiting, large 13C KIEs would be expected at both the aldehydic carbon of 5 and the α-carbon of MA, in agreement with qualitative expectations. Such large KIEs are not observed. This is compelling evidence that the aldol step is not rate-limiting, or not majorly so once the possibility of competitive rate-limiting steps is considered. McQuade’s results had already established this for the DMSO reaction.

We can now consider quantitatively whether the 13C KIEs fit with the elimination step (actually the deprotonation step of the E1cb(irr) process) being mainly rate-limiting with the aldol step being minorly rate-limiting. For the kinetic scheme of eq 2, KIE1 and KIE2 are the KIEs that would be observed if the first and second steps were purely rate-limiting, respectively. The intermediate is partitioned between a product-forming process occurring at rate rate2 and a reverse process affording the starting materials occurring at rate rate–1. It can be readily shown that the observed isotope effect KIEobs will be determined by eq 3, where the commitment factor Cf is the ratio of rate2 to rate–1. For the quantitative analysis of the methanol KIEs here, we will assume that KIE1, the isotope effect if the aldol step were rate-limiting, is approximately the B3LYP-predicted 1.042. We also assume that KIE2, the isotope effect if the elimination were fully rate-limiting, has the approximate value of 1.009 as observed for the DMSO reaction. The observed KIE would then be equal to the average experimental value of 1.0155 in unlabeled methanol when Cf = 0.245, i.e., when the second step is slower than the first step by a factor of about 4.1. Due to the H/D KIE for methanol versus d4-methanol of 3.1,[37] the Cf for the d4-methanol reaction would go to 0.08 and the expected KIEobs for the aldehydic carbon would be 1.011. This is in striking agreement with the experimental 13C KIE in d4-methanol. This agreement supports the interpretation of the KIEs in methanol as resulting from competitive aldol and elimination steps, with the latter being slower by a factor of roughly four. This rate difference defines a 0.8 kcal/mol difference between the height of the aldol and elimination barriers. Depending on the choice of assumed value of KIE1, the barrier-height difference varies by ±0.4 kcal/mol, but the observed isotope effects continue to fit well with competitive rate-limiting steps.

The Elimination Step: Eliminations in Synthesized Intermediates

To learn more about the elimination step in the MBH mechanism, we adopted the approach of independently generating an intermediate and studying its conversions under the reaction conditions. No practical synthesis of 12 itself was apparent, but the close analogue 15 was accessible by methylation of 14, the adduct of MBH product 7 and diethylamine. The salt 15 was a 2.5:1 mixture of diastereomers and was sufficiently stable to be chromatographed on silica gel. However, it could not be isolated in analytically pure form due to a slow decomposition into MA, 5, 7, and diethylmethylammonium triflate (16).

The reaction of 15 under MBH conditions using 30 mol % DABCO in d4-methanol leads to a mixture of the elimination process affording MBH product 7 plus the ammonium salt 16 and the retro-aldol process affording MBH starting materials 5 and MA along with 16. The formation of these products closely follows pseudo-first-order kinetics (see the SI). This is surprising at first glance since the formation of 16 would decrease the basicity of the solution (see the example below where this effect comes into play) but the acidity of 16 is buffered by the formation of 9, which takes up a proton. Kinetic modeling of the product concentrations versus time gave a best-fit ratio of the elimination rate constant kelim to the retro-aldol rate constant kret of 0.14:1. (This ratio was assumed to be the same for both diastereomers of 15 in order to minimize the parameters fit to experiment in the model.) From this ratio, the aldol/elimination barrier-height difference would be 1.2 kcal/mol. This result is in remarkable agreement with the 1.0 and 0.8 kcal/mol values obtained from the Eyring and 13C KIE analyses above. The three independent analyses are mutually supportive.

To approximate absolute values for kelim, the rate law was taken as rate = (kelim + kret)[15][MeO–] (see below), and the concentration of methoxide was inferred from the pKas of DABCO-H+ and 16, the methanol autoprotolysis constant,[27] and the initially measured concentrations of ammonium salts. With these assumptions, the best-fit kelim values for the major and minor diastereomers were 70 and 180 M–1 s–1, respectively. Due to potential inaccuracy in the concentration of methoxide, these kelim values have a significant potential inaccuracy, but their specific values will only be of importance in comparison with rate constants below derived from the same assumptions.

The proton-shuttle mechanism requires a free hydroxyl/alkoxide group in the pre-elimination adduct. To examine the role of the hydroxyl in the elimination, ammonium salt 18, the methoxy analogue of 15, was prepared as a 3.1:1 mixture of diastereomers by O-methylation of 7 with AgO/MeI to afford 17, addition of diethylamine, and N-methylation. Like 15, 18 could not be isolated in analytically pure form due to a slow decomposition forming 17. Based on relative chemical shifts and coupling constants (see the SI), the major diastereomer of 18 tentatively appears to correspond to the minor diastereomer of 15, and vice versa.

Under MBH conditions with DABCO in d4-methanol, 18 affords only 17 and 16. Unlike with 15, in this case there is no buffering addition reaction and the pH drops as the reaction proceeds. Accordingly, the reaction does not follow first-order kinetics, but the conversion versus time fits well with the rate law kelim [18][MeO–] (see the SI). The concentration of methoxide ion as the reaction proceeded was kinetically modeled from the concentration of ammonium salts and the same assumptions for pKas and the methanol autoprotolysis constant as used for 15. With these assumptions, the best-fit kelim values for the major and minor diastereomers were 210 M–1 s–1 and 50 M–1 s–1, respectively.

The striking result here is that the elimination occurs at nearly identical rates for the hydroxyl compound 15 and methoxy analogue 18. Reactions proceed by the fastest available mechanism, and the absence of acceleration by the hydroxyl group precludes its significant involvement in the mechanism. The combination of this observation and the absence of an H/D solvent KIE provides compelling evidence that the proton-shuttle mechanism and the oft-calculated structure 4 have no mechanistic relevance! In its place, the data support the simple acid–base protonation of 11 and deprotonation of 12.

Why? The proton-shuttle pathway might be considered to be the simplest of potential mechanisms, as it allows the direct conversion of the intermediate 11 into products with the aid of a single molecule of solvent. However, the protonic equilibrium of 11 with the solvent methanol to afford 12 involves a nearly thermoneutral[38] proton transfer between heteroatoms. Such proton transfers occur at diffusion-controlled rates, and in methanol solvent this should occur many orders of magnitude faster than the barriered abstraction of the α-C–H proton. This makes 12 an obligatory intermediate in methanol after the aldol step. The apparent simplicity of short-circuiting 12 by direct reaction of 11 is thus illusory. The elimination mechanism would proceed from 12 back through 11 only if there were some energetic advantage for such a pathway over the direct methoxide-mediated deprotonation at the α position of 12. Methanol would be recognized as a poor choice of solvent when chelation control is desired synthetically, so there is no experiential reason to expect that the chelating proton-shuttle process should be preferred. Considering this, the evidence against such a mechanism provided by the observations with 15 and 18 should not be surprising.

Proton-shuttle mechanisms have the practical advantage that they are easily explored computationally, and they have been increasingly popular observations in recent years. In contrast, two-step mechanisms involving proton transfers to and from solvent are not so readily tractable in computational studies. Such mundane acid–base mechanisms, however, have vastly greater experimental support. For example, while computational studies have routinely proposed proton-shuttle mechanisms for keto–enol and related tautomerisms,[39] experiments have strongly supported a simple acid–base mechanism.[40,41]

The Elimination Step: Transition Structures, Predicted KIEs, and Experimental Geometries

By either the proton-shuttle or acid–base pathways, the actual reaction in solution would involve an ensemble of transition states and solvation shells. No single model is likely to adequately represent either mechanism. Our computational approach to the exploration of these pathways was to obtain a variety of transition structures, a total of 27 (see the SI), by varying both the involvement of explicit methanol molecules and the involvement of the alkoxide/alcohol group of 11/12, as well as using both B3LYP and M06-2X DFT methods. Due to the unavoidably incomplete modeling of the solvation, the energies of the structures are a dubious guide to their applicability to the solution reaction. For example, the incorrect proton-shuttle mechanisms are favored over corresponding (and more accurate) simple deprotonations by 2.3–5.6 kcal/mol. In the place of a purely computational evaluation of the calculated transition structures, we use a comparison of the experimental KIEs with those predicted for the various structures. In this way, the observed KIEs can delimit some features of the experimental reaction,[42] even in the absence of a clear choice of computational model or reliable energetics.

The combination of “functional shopping” and “computational model shopping” leads to greatly varying transition structures and predicted KIEs. The range of structures and the general trends are summarized in Figure 6. All the transition structures 19 lead initially to intermediate zwitterion 20. The transition structures fall on a spectrum ranging from “early” to “late”. Deprotonation of the α-C–H bond by a “naked” methoxide anion (lacking explicit hydrogen bonding but stabilized by the PCM implicit solvent model) is relatively exothermic, and this results in early transition structures with α-C–H distances less than 1.3 Å. The predicted C-α KIEs for such structures are 1.006 or less, far from the experimental value of ∼1.013. It should be recalled that the experimental KIE is a composite arising from a mainly rate-limiting elimination step and a minorly rate-limiting aldol step. This combination was allowed for in the prediction of the isotope effects using eq 2, the predicted KIEs from Figure 5, and a Cf value of 0.245 as inferred above. We were unable to bring the composite KIE predictions into reasonable agreement with experiment with any reasonable change in either the predicted KIEs for the aldol step or the assumed Cf. This disagreement suggests that the naked methoxide/PCM implicit solvent approach is an inadequate model for the transition state.

Figure 6

Summary of transition structures for the α-C–H deprotonation. Product 20 is shown in its lowest-energy conformation. See the SI for a POV-ray structure of 20.

The methoxide ion carrying out the α-C–H deprotonation may be hydrogen bonded to additional solvent molecules in the acid–base transition structures or be hydrogen bonded to the intramolecular hydroxyl group as in the proton-shuttle transition structures. In either case the hydrogen bonding leads to later transition structures. In the M06-2X transition structures the effect is small, and a single hydrogen bond to the basic oxygen (OA in 19) leaves the transition structure relatively early (α-C–H distances of 1.29–1.30 Å). The predicted C-α KIEs for such structures, including three simple M06-2X proton shuttle transition structures, are 1.006–1.008. This disagreement with the experimental C-α KIE weighs against the accuracy of the transition structures. As an exception to this generalization, a proton-shuttle transition structure that included two methanol molecules hydrogen bonded to OB of 19 led to a predicted C-α KIE of 1.016. This is within the uncertainty of the experimental measurements. However, the predicted solvent H/D KIE for this structure was 1.59, which is inconsistent with the 0.96 ± 0.1 experimental solvent KIE.

At the opposite extreme, B3LYP calculations including two hydrogen bonds to OA led to late transition structures with α-C–H distances of 1.38–1.40 Å. The predicted C-α KIE for such structures was 1.029–1.034. This is far too high versus experiment.

In both the B3LYP and M06-2X calculations, there is a range of structures that lead to reasonably accurate predictions of the experimental KIEs. In the M06-2X calculations, structures that include two hydrogen bonds to OA lead to C-α KIEs of 1.012–1.016, C-β KIEs of 0.988–0.990, and MA carbonyl carbon KIEs of 1.009–1.012. In the B3LYP calculations, structures that include a single hydrogen bond to OA lead to somewhat less accurate but still reasonable KIE predictions: C-α KIEs of 1.015–1.017, C-β KIEs of 0.992–0.993, and MA carbonyl carbon KIEs of 1.007–1.009. The solvent H/D KIE predicted for all of these structures is in the range of 1.00–1.12; considering the general difficulties in predicting solvent KIEs, this agreement is fine.

It is disconcerting that the models leading to reasonable KIE predictions for the two DFT methods involve different levels of solvation and can involve either proton-shuttle or acid–base pathways. In this way even the limited set of transition structures that are consistent with the experimental KIEs is indecisive about aspects of the mechanism. On the other hand, there is an important commonality among the seven structures giving good KIE predictions in that they all have α-C–H distances of 1.33–1.36 Å and OA–H distances of 1.29–1.32 Å. In previous work we have shown that predicted KIEs can reflect transition state interatomic distances in a way that is independent of both the choice of theoretical method and the detailed choice of the computational model.[42] In this way, series of calculations can be used to delimit transition state distances from experimental KIEs. The α-C–H distance of 1.33–1.36 Å and OA–H distance of 1.29–1.32 Å are then a new example of the experimental measurement of a transition state geometry using KIEs. Although tunneling greatly complicates the prediction and interpretation of H/D KIEs for proton-transfer reactions, the results here interestingly suggest that 13C KIEs may be used to report on transition state distances for protons being transferred.

The Experimental Free-Energy Profile

The data discussed above along with some additional inferences define an experimental standard-state free-energy profile for the mechanism of the MBH reaction of MA with 5 in methanol catalyzed by DABCO. This profile is shown in Figure 7. It should be recognized that there are significant uncertainties associated with some of the energy values due to the nature of the estimates involved. Nonetheless, the complete free-energy profile is rooted in experimental observations.

Figure 7

Experimental and computational free energies along the MBH reaction coordinate. The black continuous line is experimental. All of the B3LYP/PCM, M06-2X/PCM, and M06-2X/SMD calculations include a full optimization with the solvent model for methanol.

In Figure 7, the basis for assigning the energies of the transition states for the aldol and elimination steps lies in the Eyring study of the reaction, the partitioning of 15, and the KIE evidence indicative of mixed rate-limiting steps. The free energy of 8 is based on the observable stability of 9 along with the kinetic acidity of 9. To get the energy of 12, we first assume that the rate constant for elimination of 12 mediated by methoxide is approximately the same as it is for the elimination in 15. This assumption seems reasonable since departure of the differing amines in the two eliminations does not occur during the rate-limiting proton abstraction. With a rate constant of ∼180 M–1 s–1 for elimination in 15, the free-energy barrier for the elimination is ∼14.4 kcal/mol. Some uncertainty in this value arises from the uncertain concentration of methoxide during the elimination reaction of 15, which in turn arises from uncertainty in the pKas of DABCO-H+ and 16 in methanol. As a check on the reasonableness of the 14.4 kcal/mol value, it may be noted that there is an identical barrier for the deprotonation of 9 by methoxide (Figure 2). This supports the assumption that the difference in the ammonium salts between 12 and 15 makes little difference in the barrier for deprotonation. If the barrier for the elimination reaction of 12 is 14.4 kcal/mol, then 12 is 6.8 kcal/mol above the starting materials.

The energy of 11 is based on an estimate of the relative pKa in water of 12 versus methanol,[38] with allowance for the differing standard-state concentrations of the two. Despite the similar standard-state free energies of 11 and 12, little 11 would be present relative to 12 because the concentration of methoxide anion would be many orders of magnitude below 1 M.

The similarity of the protonation of 20 by methanol to form 12 and the protonation of 8 by methanol to form 9 provides two semi-independent ways to estimate the energy of 20. A thermodynamic approach assumes that the acidities of 12 and 9 are equal. A kinetic approach assumes that the barrier for protonation of 20 is the same as that for protonation of 8 (13.1 kcal/mol from Figure 2). The two methods identically place the energy of 20 at 8.1 kcal/mol. The lowest-energy conformation of 20 evinces no special stabilization from either hydrogen bonding or C–OH negative hyperconjugation, but it does exhibit an extra steric interaction (depicted) not present in 8 or 12. If there is a special instability in 20, however, it does not show up in the kinetics of its formation. The barriers for formation of 8 and fragmentation of 20 are unknown.

The Computational Mechanism

Before comparing computed and experimental energetics, the choice of structures in the computational mechanism involves some subjective decisions that require discussion. In particular, we eschewed supramolecular structures including extra explicit solvent molecules. All of the intermediates and transition states along the mechanism are more strongly solvated than the starting materials and products, and the addition of explicit methanol lowers enthalpies across the board. In exploratory studies, the greatest errors, those for 9 or 12 + MeO– in any PCM calculation, were markedly decreased (by 12.2 kcal/mol with M06-2X) with explicit solvation of the MeO– by a methanol molecule. Other errors however were significantly increased, most notably those for all M06-2X free energies for 13 and 19. (The already low barrier for 19 drops by 4.7 kcal/mol with an explicit methanol molecule.) The entropy barriers for 13 and 19 are also taken substantially further from experiment by the addition of explicit solvent. The huge errors for 9 or 12 + MeO– without explicit solvation exemplify the common computational intractability of acid–base steps in mechanisms. The recognizable inadequacy of implicit solvation in these steps is not of central emphasis here, but a consequence is important. That is, due to the high computational energy of 12 + MeO–, the prior studies had no hope of identifying the correct mechanism.

The preclusion of explicit solvation excludes most variants of structure 19. The lowest-energy remaining structures are the simple proton-shuttle structures proceeding directly from 11 to 20. The latter structures were used for comparison with the experimental barrier despite their not being consistent with our experimental observations. Two arguments support the value of this comparison with an incorrect calculated mechanism. A qualitative argument is that the proton-shuttle transition states are likely to crudely resemble the actual transition state, at least to the degree that both involve an α-C–H deprotonation by a hydrogen-bonded alkoxide oxygen. A more subtle but quantitative argument is that the free energies of the proton-shuttle transition structures represent upper bounds on the free energies of the actual mechanism. The upper-bound limitation arises because reactions must occur by the lowest-energy mechanism so any incorrect mechanism in reality must be higher in energy.[43]

Computational Energetics. The Entropy Problem and an Approximate Solution

The free energies to be presented were calculated in two ways. In the first way (eq 4), the free energy Gtot was derived by adjusting the raw harmonic entropy (Sharm,1 atm) to a 1 M standard state (except for the neat standard state of methanol) and correcting for entropies of mixing (Smix, which allows for enantiomeric and other low-energy conformations, see the SI) and symmetry (Ssym). These well-known corrections are significant here, though they are almost universally neglected in the literature. In the second way (eq 5), the free energy G50% was calculated after halving the raw entropy.

The calculation of entropy changes in solution has been considered to be a substantial problem since differential solvation entropy may play a large role. The purpose of considering G50% is that it explores an example of literature methods that have been used to correct for the perceived problems in the calculation of entropies in solution. Such problems are illustrated by a comparison of the experimental entropy of reaction for formation of 9, −24.9 e.u., with calculated harmonic entropies of reaction after adjustment to a 1 M standard state ranging from −34.6 to −41.0 e.u. The calculations thus overestimate the decrease in entropy by approximately 10–16 e.u. The errors appear larger for the third-order aldol and fourth-order elimination steps, overestimating the entropic barriers by 46–56 and 25–34 e.u., respectively. (The frequent mistake of failing to convert to a 1 M standard state would aggravate these errors.) Such overestimates of the entropy loss in association reactions in solution have been a common observation. Errors in calculated entropies in solution could have many sources, but they have most often been attributed vaguely to restrictions on translational and rotational degrees of freedom.

Rigorous approaches to the calculation of solution entropy changes can be difficult to apply in general,[44] so diverse estimation tactics have arisen for dealing with the entropy problem. Some computational papers reduce the entropy by 14.3 e.u. (this is R ln(1354), which is the ratio of 55 M to 1 atm).[45,46] Some use related adjustments differing by solvent.[47] Some cut the entropy by 50%.[48] A similar correction was employed by Aggarwal and Harvey for the MBH reaction.[11] Some employ more elaborate entropy adjustments based on ideas originated by Wertz.[49,50] Some argue that translational and rotational entropy in solution should be ignored altogether.[51] Considering arguments in the literature over the nature of the entropy error[52] and the absence of a general theoretical basis that would apply to diverse solvents and solutes, the physical basis for assuming any particular entropy error is unclear. Notably, some reactions exhibit greatly decreased experimental entropy barriers in solution versus the gas phase[52] while others exhibit no decrease at all.[53] The range of entropy corrections commonly used in the literature provides roughly a 100,000-fold range of choices for predicted equilibrium or rate constants for bimolecular association reactions, and a 1010-fold range of choices for trimolecular association reactions. The 50% entropy correction explored here is typical and illustrative of the range of defended computational predictions.

There is, however, a fundamental problem with the usual consideration of entropy in computational studies of solution reactions. Real free energies of solvation are temperature dependent. Experimental entropies are determined from the temperature dependence of rates or equilibria, so they reflect the temperature dependence of solvation free energies. Neither gas-phase calculations nor most continuum solvation models, e.g. PCM and SMD, allow for the temperature dependence of the solvation free energy in any way. As a result, the entropies that arise should not be expected to be comparable to experimental values. This is equally true for experimental versus calculated enthalpies.

In an attempt to address this problem, we made use of the temperature-dependent SM8T solvent model of Cramer and Truhlar. The SM8T-calculatd free energies of solvation (ΔGsol) at two temperatures (here T = 308.15 and T0 = 298.15) were used to define a phenomenological entropy of solvation (ΔSsol, which also incorporates a term for the heat capacity of solvation) as in eq 6. The substantial ΔSsol values, e.g., −69 e.u. for 12, were used in combination with Smix and Ssym to correct the total entropy Stot for each molecule as in eq 7. The enthalpy Htot is then calculated from Gtot and Stot or it may be viewed as the result of a counterbalancing adjustment to Hharm (eq 8). This simple but previously unused process has the effect of partitioning the free energy derived from PCM or SMD calculations into entropy and enthalpy components suitable for comparison with experiment. We would emphasize that this partitioning does not change Gtot.

For the four cases with experimental entropies (7, 9, 13, 19), the uncorrected calculated entropy at a 1 M standard state is consistently far too low, with a mean absolute error (MAE) across all computational methods of 25 e.u. The combination of the 50% entropy correction, Smix, and Ssym reduces the MAE to 18 e.u. This literature adjustment thus improves predictions but the errors are still large, including an average error of 41 e.u. (9 orders of magnitude) for 19. The combination of ΔSsol, Smix, and Ssym incorporated in Stot reduces the MAE to only 9 e.u. Moreover, 9 of the 12 total entropy predictions err by only 2–6 e.u., which is similar to or less than the experimental uncertainties. Only the entropy of 13 remains subject to a substantial error of 20–30 e.u. in predicted entropies of activation, and the experimental entropy for this partially obscured aldol step is the most uncertain of the experimental values.

An important observation is that after this best-practical entropy calculation, the errors in ΔHtot far exceed those from TΔStot, even for 13. (See Figures 2 and 7 for the individual comparisons.) As far as can be discerned from the checkable cases, the large errors in the free energies to be discussed have their origin in the miscalculation of enthalpy, not entropy.

Comparisons of Experimental and Calculated Energetics

Figures 2 and 7 compare the experimental free energies associated with the MBH mechanism with free energies derived from the B3LYP/PCM, M06-2X/PCM, and M06-2X/SMD calculations including the full harmonic entropy (ΔGtot), free energies from B3LYP/PCM and M06-2X/PCM calculations including a 50% entropy correction (ΔG50%), and some additional free energies based on high-level calculations. The most striking feature of the figures is the sheer range of computational predictions. The predicted equilibrium constants for the simple formations of ammonium salt 9 and neutral 7 vary by 12–14 orders of magnitude. The predicted equilibrium constant for formation of 11 varies by 26 orders of magnitude. The predicted energy of transition state 19 varies by 48 kcal/mol, equivalent to a range that is 35 orders of magnitude. The range of predictions would be even larger if free energies from M06-2X/SMD calculations with a 50% entropy correction (see the SI) were included in the figures. These last predictions might be defensible on a literature basis but they are too poor to warrant presentation, outside of the interesting if appalling prediction that intermediate 20 would be more stable than product 7.

The figures exhibit most obviously the exceptional error of B3LYP ΔGtot calculations and the M06-2X ΔG50% calculations. Considering the known problems of B3LYP with the energies of σ bonds relative to π bonds,[54] the B3LYP energetics might be considered to be a straw man, despite their use in many of the MBH mechanistic studies. However, the large errors in the B3LYP energies are not solely the result of an expected σ/π energy error. For the proton-shuttle transition state found by Aggarwal and Harvey, the G3MP2 correction to the energy was only 8.5 kcal/mol.[11] Applying this correction to the 44.9 kcal/mol B3LYP barrier would leave it 15.2 kcal/mol too high, still underestimating the rate by 11 orders of magnitude.

The spectacular failure of M06-2X ΔG50% would seem less easily anticipated. The M06-2X method would be expected to provide fairly accurate energetics. The literature support for the use of an entropy correction is widespread, and a similar entropy correction was employed by Aggarwal and Harvey. In the absence of experimental observations for the reaction at hand, there would be no compelling reason to discount the computational results. However, the M06-2X/50% entropy calculations err by 24.2 kcal/mol, 18 orders of magnitude, on the rate of the elimination step. This leads to the warped prediction that the formation of 8 is the rate-limiting step. This is clearly providing no useful information.

The errors in the M06-2X/PCM ΔGtot and B3LYP ΔG50% calculations are still large but are often much less extreme. If we ignore the huge errors for 9 or 12 + MeO– without explicit solvation, the largest remaining error is seen for 11 (5 and 11 orders of magnitude for the M06-2X/PCM ΔGtot and B3LYP ΔG50 energies, respectively). However, the errors are notably inconsistent. The B3LYP ΔG50% calculations have errors that range from +2.2 to +15.1 kcal/mol, while the M06-2X ΔGtot calculations have errors that range from +6.7 to −4.2 kcal/mol (ignoring 9 or 12 + MeO–). This inconsistency is vexing. On a detailed level, it leads to a number of strange predictions. One example is that 8 and 20 are structurally rather similar, yet the M06-2X calculations err in their relative free energies by 6.1 kcal/mol. Another example is seen for the aldol reaction of 8 with 5. With an M06-2X-predicted ΔG⧧tot = 5.1 kcal/mol, this step would be expected to occur at a nearly diffusion-controlled rate when its actual rate constant is only 104 M–1 s–1. In general, inconsistent errors can readily lead to qualitative errors in the prediction of a mechanism. In the MBH reaction, the M06-2X ΔGtot calculations do not predict the correct rate-limiting step, and the B3LYP ΔG50% calculations do not predict that the reaction would succeed at all.

The inconsistency of errors also precludes any inference of accuracy for one calculation based on that for another. The M06-2X ΔGtot calculations perform well for the equilibrium constant for formation of 7, but they are off by 4.1 kcal/mol in the formation of 9. The B3LYP ΔG50% calculations are strikingly good for the equilibrium formation of 9 (due to a fortuitous cancellation of entropy and enthalpy errors) but are poor for the formation of 7.

High-Level Energy Corrections

In large systems such as the MBH reaction, a potential computational approach to the minimization of errors involves the comparison of the applicable DFT methods to high-level gas-phase calculations for a model chemical system. The DFT errors in the model system are then used to correct the calculated energetics for the full system. We consider here the impact of this energy-correction process. It should be recognized that the high-level correction process is not readily applicable to most of the MBH mechanism. For example, structures 8, 11, and 20 are zwitterions for which simple and reasonably analogous gas-phase models are unstable.[55]

To address the expected bias in the B3LYP energetics against the formation of 9, we examined the model addition of NH4+ to MA to afford 21. In gas-phase calculations, the composite G3B3 method placed the formation of 21 as 1.5 kcal/mol more exothermic than B3LYP did. Applying a 1.5 kcal/mol correction factor to the B3LYP ΔGtot results with MA/DABCO-H+ versus 9 decreases the free energy error modestly (Figure 2). Surprisingly, the G3B3 correction provides a much greater improvement to the M06-2X ΔGtot prediction for 9, bringing the error down to only 0.7 kcal/mol.

The error in the B3LYP ΔGtot for the formation of product 7 from 5 and MA was particularly striking since it leads to a miscalculation of an equilibrium constant between ordinary neutral molecules of over 9 orders of magnitude, and because this error can be attributed entirely to ΔHtot, not ΔStot. In an attempt to correct for this error, the formation of 7 was modeled by the reaction of ethylene with acetaldehyde to afford 3-buten-2-ol (22). The G3B3 correction of 5.2 kcal/mol significantly improves the B3LYP ΔGtot predicted energetics, but the equilibrium constant would still be off by over 5 orders of magnitude. Interestingly, in this case the G3B3 correction makes the M06-2X ΔGtot free energy worse by 2.2 kcal/mol.

The very largest of the errors in the full system occur for the formation of 12 + MeO–. It was of interest whether corrections based on high-level calculations could decrease these errors. The reaction was modeled by the reaction of MA with ammonia and formaldehyde to afford 23. The correction lowers the error in the B3LYP ΔGtot modestly, from 39.7 to 36.3 kcal/mol. However, the correction raises the error in the M06-2X/PCM ΔGtot for 12 + MeO– by 7.5 kcal/mol. As a result, the equilibrium constant for formation of 12 would be miscalculated by 22–29 orders of magnitude after G3B3 correction. Clearly, the problem here lies with the inadequacy of the solvent model. The high-level correction provides only false reassurance.

General Perspective

The MBH reaction in methanol is in some ways a challenging system for computational mechanistic chemistry. Each intermediate is either zwitterionic or charged, the effect of the polar solvent in stabilizing the charges is large, and the error inherent to any implicit solvent model in estimating the charge stabilization by solvent would be expected to be substantial. In addition, the second-order and third-order nature of most of the intermediates and transition states along with the minimally fourth-order nature of the rate-limiting step inflates the role of entropy in the relative free energies along the mechanistic pathway, and this maximizes the potential error due to the misreckoning of entropy. Error in the MBH mechanism should not be a good exemplar for error in some simpler reactions, such as nonpolar unimolecular pericyclic rearrangements. On the other hand, all of the intermediates and transition states are closed shell species without low-lying excited states and there is no reason to expect that electronic structure methods should intrinsically lead to sizable errors in the MBH reaction. Overall, the challenges imposed by the MBH mechanism would not appear to be greater than those seen in a large portion of computational studies of reactions in solution.[56]

There is no practical limit to the number of alternative computational method, basis set, entropy calculation, and solvent model combinations that could be applied to the MBH mechanism. It is inevitable that some subset of the possible computational approaches will provide a more accurate prediction of the free energy surface for the MBH reaction. Such accuracy is meaningless unless the computational approach would reliably and foreseeably make accurate predictions. While any popular computational approach has its virtues, none are used more often in the recent literature than those applied here. Few would provide a cogent reason to expect more accurate results for the MBH reaction. The M06-2X results using the SMD solvent model were originally explored as a logical approach to decreasing error. The error in the energy of 12 + MeO– decreased greatly relative to the PCM results. Other errors increased greatly, particularly those for 13 and 20.

Additionally, each structure in Figure 7 was reoptimized in PCM calculations with the larger 6-311+G** basis set. The mean absolute deviation from experiment increased by 1.8 and 0.2 kcal/mol for the B3LYP and M06-2X functionals, respectively.

Finally, we explored the full MBH mechanism using a series of computational methods including empirical dispersion corrections,[57] including B3LYP-GD3/6-31+G**/PCM, B3LYP-GD3/6-31+G**/SMD, M06-2X-GD3/6-31+G**/SMD, PBE0-GD3/6-31+G**/SMD, and wB97xD/6-31+G**/SMD methods. The results from these calculations are given in the SI. The notable success of inclusion of the D3 dispersion correction is that the greatest errors of the B3LYP ΔGtot calculations are greatly curbed. For example, the barrier associated with 19 is reduced from 44.9 kcal/mol in the B3LYP/PCM ΔGtot calculations without a dispersion correction to 24.0 kcal/mol including the dispersion correction. It might be argued in general that the combination of dispersion corrections and the SMD solvent model would minimize error, and this combination in fact leads to further improvement in the B3LYP energies. However, the same combination uniformly increases the errors in the M06-2X ΔGtot results. Overall, ignoring the pernicious case of 12 + MeO–, the uncorrected M06-2X/PCM results in Figure 7 were the most accurate. Four general observations should be noted. First, every method explored led to at least one free-energy error of more than seven kcal/mol, or 5 orders of magnitude. Second, every method erred on the relative energies of 11 and 12 by more than eight kcal/mol. Third, no method was most accurate in its prediction of more than one of the experimental energies. Fourth, every set of ΔG50% free energies included one or more errors in excess of 15 kcal/mol.

Computational methods are simply scientific models. Any model makes some inaccurate predictions but models can retain utility despite significant propensities for inaccuracy. Inaccurate predictions aid the choice of models for future predictions. Because of this, the central scientific problem in the computational study of the MBH mechanism is not the inaccuracy of the predictions. Rather, it is the absence of any particular prediction at all. Fully defined computational methods (including the choice of basis set, entropy calculation, and solvent model) of course make quite specific predictions. However, there is neither a consensus best-choice method nor a common view on the right way to choose a method. When evaluable, the most accurate choice varies with the system at hand. In the MBH reaction, defensible and expectantly publishable choices of computational approaches lead to predictions of the facility of the proton-shuttle process that vary by 35 orders of magnitude in the stability of 19, while also diverging in the geometry and preferred stereochemistry of transition state 13. This variance is in practical terms indistinguishable from making no prediction. In addition, studies of the MBH mechanism have not been considered falsified by extreme inaccuracies in predictions.[13b,14b,15b] In the terminology of Pauli, computational mechanistic chemistry is “not even wrong” about the MBH mechanism.

A less bleak view of the utility of computations in the study of the MBH mechanism can be built around the argument that the experimental observations for a reaction can and should be used in the choice of theoretical methods used in the study of that reaction. Even a single comparison with experiment is very helpful; the poor predictions of the overall barrier for the MBH reaction (which is easily estimated without any detailed kinetic study) would allow one to exclude the otherwise defensible B3LYP ΔGtot and M06-2X ΔG50% calculations. If one simply requires consistency with a second experimental observation, i.e., that the reaction of MA with 5 proceeds, then the M06-2X ΔGtot calculations would be chosen. This is a tremendous advance over the incredible range of predictions that might be obtained in the absence of consideration of experimental observations, in part because the M06-2X ΔGtot calculations provide the best overall prediction of the free-energy profile but more importantly because the delineation of a specific method leads to specific and testable predictions. However, the M06-2X ΔGtot calculations qualitatively mispredict the rate-limiting step, quantitatively mispredict the relative energy of 8 versus 20, and substantially underpredict the enthalpies of the key transition states so this process certainly does not preclude incorrect predictions.

Conclusions

In outline, the very shortest mechanism for the MBH reaction here would involve only four steps through 8, 11, and 20. However, the uncatalyzed conversion of 11 to 20 is not tenable owing to the well-known slowness of direct 1,3-intramolecular proton transfers.[41] Allowing for this, the simplest viable “electron-pushing” mechanism would proceed in five steps through 8, 11, 12, and 20, with the solvent methanol accelerating the reaction by mediating the acid–base steps. It is exactly this mundane mechanism that is supported by the observations here.

A key conclusion is that the conversion of 11 to 12 occurs by acid–base chemistry and not by a proton-shuttle process. The absence of a proton-shuttle mechanism is supported by a 0.96 ± 0.1 solvent H/D KIE, and it is strongly supported by the nearly identical rates for eliminations of 15 and 18. This conclusion is in contrast to the seven computational mechanistic studies that had previously considered this issue. The preferred consideration of computationally tractable mechanisms over less tractable alternatives is common, and the results here underscore that there is no scientific basis for this preference.

Our experimental observations define a nearly complete free-energy profile for the mechanism. A series of observations support the involvement of competitive rate-limiting steps, and this allows both barriers to be determined. The proton-transfer step forming 20 is the primary rate-limiting step at 25 °C, but the aldol step forming 11 is partially rate-limiting, and it becomes the primary rate-limiting step at low temperatures. Other observations define the energy of 9 and 12, indirectly delimiting the energies of 8, 11, and 20.

The general outline of the MBH mechanism as an addition/aldol/elimination sequence was understood from experimental observations before any computational mechanistic studies. The McQuade mechanism employing a second molecule of aldehyde to facilitate the elimination arose from experimental studies. The suggestion by Aggarwal and Lloyd-Jones that hydroxylic compounds accelerated MBH reactions in a similar way, though not supported here in its details, arose from experimental observations. Though never emphasized, computational studies successfully recognized the E1cb nature of the elimination step. This conclusion however would have been clearly anticipated from experimental studies.[58] The more primordial currency of information provided by computational studies consists of the geometries and energies of intermediates and transition states along the mechanism. Except for the ultimately irrelevant proton-shuttle and direct 1,3-proton-transfer transition states, the various geometries received little discussion in the published studies. The results here with 19 and 13a–c highlight the large variations in geometries and changes in the preferred diastereomer that occur with changes in the theoretical model and computational method. We have discussed in detail the problems with the computed energies. Overall, it is not clear to us that any significant accurate information that was not already apparent from experiment either has been, or could have been, reliably garnered purely from computations. In the absence of any consideration of experimental observations in the MBH reaction, defensible computational studies could have made an exceptional diversity of predictions, many of which would have been absurd. In the actual MBH case where much was known experimentally, computational predictions that were consistent with experiment were emphasized while those inconsistent with experiment, such as the B3LYP findings of astronomical barriers and that the product was less stable than the reactants, were ignored. The computational studies then highlighted one essentially pure prediction—that of the proton-shuttle process—and that prediction was incorrect.

From a more positive perspective, our results lead to some specific recommendations for computational mechanistic studies. The general and most often spectacular failure of the ΔG50% energies suggests that arbitrary entropy “corrections” should be abandoned. The interpretation of published studied including these corrections should be approached with great care. The failure of the proton-shuttle mechanism for the MBH mechanism suggests that such pathways should always be weighed carefully against simple acid–base mechanisms. The consistently large energetic error in the proton transfer converting 11 to 12 should be noted, as well as a similarly large error in the analogous proton transfer between alcohols in a previous study.[11] The experimental literature on acidity and proton-transfer rates is massive, and we would suggest that the facility of proton transfers is often best evaluated on the basis of this literature instead of direct calculation. Finally, the errors in relative energetics seen here should be considered in the credence given to the assignment of mechanisms and rate-limiting steps from computational mechanistic studies.

Computations aided significantly in the mechanistic interpretation of the experimental 13C KIEs in terms of a commitment factor and the mixture of rate-limiting steps involved. Computations intriguingly also provide a detailed, model-independent geometrical interpretation of the 13C KIEs in terms of interatomic distances in the elimination transition state. Regardless of the associated uncertainty, computations remain the only available handle on the transition states for formation of 8 and fragmentation of 20. Overall, the combination of experimental and computational studies provides a full mechanistic pathway for the MBH reaction including details that would be impossible to discern from either alone.

The scientific approach taken here has been that of a case study, and as such it suffers from the general limitations of case studies. The most important of these is the problem of generalization of the results to a broader swath of cases. The problems in the computational study of mechanisms encountered in the MBH reaction certainly cannot be used to paint all computational mechanistic studies. Many, either by simplicity or carefully designed use of the computations, would not be susceptible to the difficulties encountered here. At least, however, it would seem that studies of complex multimolecular polar reactions in solution should be undertaken and interpreted only with extreme care. The strength of a case study is that it identifies problems for consideration in other cases, and the results here suggest a variety of issues that should be carefully considered in the execution and interpretation of computational mechanistic studies.