Article pubs.acs.org/Organometallics
Utility of the Nudged Elastic Band Method in Identifying the Minimum Energy Path of an Elementary Organometallic Reaction Step Kate E. McPherson, Libero J. Bartolotti, Andrew T. Morehead, Jr.,* and Andrew L. Sargent* Department of Chemistry, East Carolina University, Greenville, North Carolina 27858, United States S Supporting Information *
ABSTRACT: The nudged elastic band (NEB) method has been used to re-examine the oxidative addition step in the classic rhodium-catalyzed hydroacylation reaction. Numerous additional intermediates were found on a pathway that is lower in energy than that originally reported. This study illustrates the potential limitations of the traditional approach to locating transition states in which a chemical-intuition-guided linear or quadratic synchronous transit is used to estimate the high point of the reaction trajectory. Utilization of that approach constrains the search so that the actual transition state and other possible intermediates and transition states may be missed.
■
INTRODUCTION Transition-metal-catalyzed reactions are widely utilized in organic synthesis,1 yet the mechanisms of these reactions are notoriously difficult to examine experimentally, as the concentrations and kinetic lifetimes of the intermediates are often minimal. Computational modeling can lend insight, but can also be very difficult, as multiple competing pathways may have similar energetic profiles.2 The sheer number of transition states that reside on the potential energy surface of such reactions is staggering, and navigating a path across the surface that agrees with experimental data involves a great deal of decision making. In this paper, we review the standard decisionmaking process that was involved in our previously reported3a mechanistic analysis of the Bosnich rhodium-catalyzed hydroacylation (HA) reaction4 and compare it to a new approach that utilizes the nudged elastic band (NEB) method5 for charting a path across the potential energy surface. Our previous model of the overall reaction pathway for the rhodium-catalyzed intramolecular HA of 4-pentenal is illustrated in Scheme 1. As proposed by Bosnich, productive hydroacylation was found to proceed from intermediate 1 via oxidative addition to form 2a, followed by alkene insertion, resulting in intermediate 3, and reductive elimination to produce the coordinated ketone product 4. The competitive pathway 1 → 2b → 5 → 6 → 7 → 8 resulted in an inactive complex and catalyst death. The relative energetics of the entire process was consistent with experimentally determined relative rates and isotopic labeling studies. Important observations included (a) the stereochemistry of the oxidative addition product determined the pathway followed, (b) coordination of a second Lewis base to 3 played a critical role in the partitioning between hydroacylation and decarbonylation, and (c) reductive elimination was rate-limiting. © XXXX American Chemical Society
Many organometallic transformations masquerade as deceptively simple processes. Consider the “fundamental” step for the oxidative addition of the C−H bond in this reaction. Although there are multiple pathways to arrive at the goal, the simplest involves the direct cleavage of the C−H bond. Sigma and pi complexes are often postulated as intermediates in C−H activation reactions,6 but are often difficult to locate in computational models of the minimum energy path (MEP) for such reactions.7 In our earlier work, the transition states for the direct C−H cleavage pathway connecting 1 and 2a or 2b were particularly difficult to identify. We used a traditional approach that employed a linear or quadratic synchronous transit based method to estimate the high point of the reaction trajectory.8 The selection of this point was tempered by our chemical intuition that suggested, in this case, the peak should coincide with the point of C−H bond cleavage. A local surfacewalking algorithm was utilized to locate the transition state, and considerable computer time was required to traverse what seemed to be numerous soft modes in the surface and arrive at a saddle point. Satisfying to us at the time, the saddle point appeared to coincide with C−H cleavage, and the subsequent frequency calculation and animation of the imaginary frequency appeared to verify the connection between this point and the adjacent minima. One of the key (and very common) assumptions we made was that the lowest energy pathway would proceed through the lowest energy intermediates. As can be seen in Scheme 1, there are four possible diastereomeric intermediates (2a, 2b, 2c, and 2d) formed by oxidative addition. Following the usual practice, we chose the two lowest. While we were aware that a sigma or Received: March 22, 2016
A
DOI: 10.1021/acs.organomet.6b00236 Organometallics XXXX, XXX, XXX−XXX
Article
Organometallics Scheme 1. Morehead and Sargent Minimum Energy Pathways for Hydroacylation and Decarbonylationa
a
Energies shown are ΔE relative to structure 1 = 0.0 kcal/mol. relativistic corrections,3f etc.; see Supporting Information for additional details).
pi complex could be present along the pathways, the structures we identified were higher in energy than the lowest oxidative addition product and were therefore ignored. In this paper, we re-examine the oxidative addition step of the HA reaction, which, with its multiple competing pathways, is particularly illustrative of the challenges facing computational modeling studies, and we use it herein as an example of the utility of NEB in mapping out a plausible reaction mechanism. Once two stationary points have been located by optimization, NEB finds an MEP that connects the two minima. NEB is a chain-of-states procedure in which atoms are mapped in a oneto-one method between the starting point and the ending point and a series of interpolated geometries (called images) are generated that interconvert the two species. The total forces are evaluated for all images (except the end points) and subsequently projected into two forms that are applied to the images: the so-called parallel forces that maintain an even distribution of the images along the MEP and the perpendicular forces that are optimized to a minimum. A common analogy of the method is that of throwing a rope over a mountain pass and dragging it back and forth to locate the lowest elevation path. The perpendicular forces will pull the images toward any minima present along that path, and therefore the NEB method significantly improves the odds of locating additional stationary points that might otherwise have been overlooked. An added bonus is that the algorithm guarantees that the transition state is connected to the adjacent local minima.
■
■
RESULTS AND DISCUSSION Motivated by our recent success using the NEB method to analyze reaction mechanisms of related systems9 (manuscripts in preparation) the Bosnich HA oxidative addition step was reexamined. NEB minimum energy paths to all four oxidative addition intermediates (1 → 2a, 2b, 2c, and 2d) are shown in Figure 1. Along the path shown in red to 2a several important features are visible, most critically that the ∼27 kcal/mol barrier
METHODS
Figure 1. Overlaid NEB pathways for 1 → 2a, 1 → 2b, 1 → 2c, and 1 → 2d. Note that the x-axis is the 17 movable images that are evenly spaced between two terminal intermediates and that being in the same position on the x-axis does not imply that the two images are equivalent.
Although NEB can be interfaced with any electronic structure program, the results presented in this report utilize the same computational methods used in our previous report3a (fully numerical DFT program DMol3,3b,c BOP functional,3d,e DNP basis, scalar B
DOI: 10.1021/acs.organomet.6b00236 Organometallics XXXX, XXX, XXX−XXX
Article
Organometallics is much higher than the previously located transition state of 19.7 kcal/mol. While the C−H cleavage step occurs at the highest point of the MEP, it is not the same structure as was found in the original report. Additionally, there is an apparent discontinuity manifested by a large jump between adjacent images at approximately point A. Experience suggests that this is the result of two pathways crossing. Furthermore, there is a plateau (point B) that suggests the presence of another potential intermediate, which, upon further investigation, turned out to be a pi complex of the carbonyl. It is noteworthy that the originally identified transition state is not actually on this or any of the other MEPs shown in Figure 1. The NEB reveals that the path directly connecting 1 → 2a is too high to be consistent with the experimental evidence that reductive elimination is the rate-limiting step. Surprisingly, the lowest energy path, according to NEB, is through intermediate 2c, which was originally not considered as part of the productive hydroacylation pathway, as the local minimum initially located was 8.0 kcal/mol higher in energy than 2a. A slight difference between minimum 2c (+15.1 kcal from 1) and the one reported in the original study (+13.2 kcal) is due to conformation: the former adopts local C2 symmetry of the phenyl groups, while the latter has Cs symmetry. The NEB for the 1 → 2c process is shown in blue in Figure 1 and possesses several interesting features. As before, the plateau at point C indicates the presence of a pi complex (2π), which then rearranges to a sigma complex 2σ (point D). (Note that the blue/green color combination for structures 2π and 2σ in subsequent figures and schemes reflects that both blue and green paths from Figure 1 share these common intermediates; see discussion below of the path to intermediate 2b. Although the structure of the 2σ complex is η1 H−Rh rather than the more traditional side-on η2 CH−Rh, the two electrons that form this three-center two-electron bond come from the C−H sigma bonding pair; therefore, we shall refer to this bonding mode as a sigma complex for the purposes of this discussion.) The sigma complex 2σ then undergoes oxidative addition to generate 2c. While not shown in Figure 1, application of similar methodology to the 2c → 2a rearrangement reveals a more complex process than expected, as 2c does not directly give 2a. First, there is inversion of the aldehyde backbone of the substrate ligand (2c → 2c′; see Figure 2). Next is a pseudorotation of the substrate to move the hydride and acyl groups from the apical and equatorial positions, respectively, of the square pyramid to the equatorial and apical positions, respectively (2c′ → 2a′). Lastly, the phosphine phenyl groups rotate (2a′ → 2a) to give the structure that goes on to hydride insertion. Performing these steps out of sequence or in a concerted manner corresponds to higher energy pathways. The oxidative addition leading to intermediate 2b and, ultimately, decarbonylation and catalyst death goes through a similar process initially. From the same sigma complex, 2σ shown on the green path in Figure 1 as point D (which comes from the pi complex C), the C−H bond rotates in the opposite direction relative to the rhodium center but does not put the hydride apical upon cleavage. If it did, we would end up with high-energy minimum 2d, whose structure was not located in our original report. Instead, the cleavage step results in an apical acyl group corresponding to structure 2b′ (Figure 2). Relatively minor rearrangements (phosphine phenyl group, aldehyde backbone rearrangement) generate 2b. As was the case along the pathway leading to 2a, the TS structure on the path to 2b found by traditional methods does not actually exist on the
Figure 2. Intermediates and transition-state structures for the competing 1 → 2a and 1 → 2b oxidative addition pathways highlighting the key geometric change between each step.
MEP found by NEB. The new activation barrier for this step is lower in energy and is guaranteed by the NEB algorithm to be connected to the adjacent minima. Optimized energies for both pathways are shown in Scheme 2, along with the energies of the transition states located during the initial study for comparison.10 Several of the other steps along the HA and decarbonylation pathways were examined with the NEB method in our initial report, and those that were not were investigated anew herein. In the new analysis, the ratedetermining step for the two complete reaction pathways that include intermediates 2a and 2b remains the reductive elimination and the carbonyl deinsertion, respectively, and thus the conclusions of the earlier publication remain consistent with the experimental data. Still, the new features of the reaction pathways change our understanding of the MEP for each of the oxidative addition pathways. These changes and their implications are discussed in further detail below. C
DOI: 10.1021/acs.organomet.6b00236 Organometallics XXXX, XXX, XXX−XXX
Article
Organometallics
Scheme 2. MEPs for Oxidative Addition Leading to 2a and 2b and Comparison with the Transition States Originally Located via Synchronous Transitsa
a
Original TS for 2a shown in red and for 2b in purple. ΔE values are listed in kcal/mol.
that higher energy intermediates are avoided along the lowest energy reaction path. Note, however, that the lowest energy nonproductive (decarbonylation) path does not include intermediate 2d but evolves, instead, through structure 2b′. This stands in contrast to the productive pathway, in which the direct C−H cleavage from 2σ to structure 2a′ or 2a necessitates the simultaneous inversion of the ethylene backbone and rotation of the alkene while cleaving the C−H bond (i.e., the path shown in red in Figure 1), making the path through 2c favorable in comparison. As in the productive pathway, the decarbonylation route has several low-energy backbone/ligand rearrangements that are necessary before alkene insertion can occur. For 2b′, rotation of the phenyls around the P−C bonds to give 2b″ is necessary for the next step, and inversion of the aldehyde backbone is required to give 2b. Interestingly, and like what we described for structure 2c above, intermediate 2b (+8.3 kcal/mol, Scheme 2) is not the same as the one originally reported (+7.8 kcal/ mol); the local minima differ by minor phenyl group rotations. We found that the new, slightly higher energy structure 2b can proceed to hydride insertion without rearranging to the lower energy isomer.
The presence of pi complex 2π and sigma complex 2σ plays a critical role in lowering the energy of the C−H cleavage. Intermediate 2σ (the C−H sigma complex) undergoes C−H bond cleavage with a barrier of 19.1 kcal/mol and results in the formation of 2c. It is quite possible that the actual process avoids the metastable 2σ species and undergoes oxidative addition directly from 2π, but the stabilization provided by 2σ is apparent in the relative rapidity of the C−H activation process. As mentioned previously, 2c cannot pseudorotate to 2a readily, but inversion of the aldehyde backbone to 2c′ proceeds with a low barrier; then 2c′ pseudorotates to produce 2a′. Phosphine phenyl rotations interconvert 2a′ to 2a. While relatively low energy processes, such small motions may have significant implications in catalyst design and/or substrate specificity. Note that the formation of the pi complex, formation of the C−H sigma complex, and the C−H cleavage are all separate steps. Movement that combines multiple steps into single concerted motions, while appealing, tends to result in higher barriers. From these data we conclude that nature appears to favor multiple, simple, elementary steps along a reaction pathway, rather than fewer, complex, concerted motions. The nonproductive pathway leading to catalyst death initially proceeds through several of the same steps. The branch point is 2σ, where rotation of the aldehyde C−H bond in opposite directions relative to the bound alkene can generate 2c (H apical, acyl equatorial) or 2b′ (H equatorial, acyl apical). The very high trans influence of the acyl group (even larger than that of hydride) means that structures 2c and 2d, with the equatorial acyl, are substantially higher energy (15.1 and 17.8 kcal/mol above baseline, respectively) than the corresponding 2a and 2b (5.2 and 8.3 kcal/mol), where the acyl is apical in the five-coordinate square pyramidal structures. In this context, it is interesting that the productive pathway therefore includes intermediate 2c. The result highlights the frequently held myth
■
CONCLUSIONS The nudged elastic band method was used to re-examine the oxidative addition step in the Bosnich HA reaction. Productive (leading to product) and nonproductive (to decarbonylation) pathways were found that were lower in energy than those identified in our initial report, and multiple additional intermediates were discovered along both routes. While the new information garnered from the NEB analysis did not alter the qualitative conclusions from our initial report, the importance of the study lies in how it illustrates certain limitations of exploring a complex reaction surface while depending on chemical intuition and traditional transition-state D
DOI: 10.1021/acs.organomet.6b00236 Organometallics XXXX, XXX, XXX−XXX
Organometallics
■
ACKNOWLEDGMENTS East Carolina University and the National Science Foundation (CHE-0415484, A.T.M.; CNS-0619285, A.L.S.) supported this work.
location techniques. For example, our analysis of a related HA system9m (manuscript in preparation) found that the highenergy point of the OA step does not coincide with C−H cleavage. Traditional methods of locating the TS structure in this case misrepresent the mechanistic step to the extent that the qualitative conclusions based on the respective approaches differ. In this context, we were fortunate that the results from our initial report correlated with experiment. Additionally, since the steric demands of the pi and sigma complexes identified herein are significantly different from those of the standard Rhcarbonyl complex 1, and the rotation of P−Ph bonds and substrate backbone feature prominently in the mechanism, these motions must be considered in any design changes of the catalyst. It is important to remember that in multistep reactions, such as the one examined herein, OA and RE, while not microscopic reverses of each other, still have many similar steric and electronic constraints. Any changes in either one is likely to have a cascade effect on all the other steps of the catalytic cycle and will thereby significantly impact the entire reaction. That the NEB method was helpful in elucidating previously unknown fine structure in the OA mechanism highlights its utility in catalyst design. Lastly, the presence of additional intermediates impacts the energetics of the C−H activation step and therefore influences the modeled lifetimes of the intermediates, the composition of the associated equilibria, and the rates of the accompanying elementary steps. The NEB method has been shown to be a valuable tool for mapping reaction pathways and provides an alternative approach to locating transition states. The algorithm facilitates the identification of additional intermediates that may be present along the reaction pathway; viewed in a different light, it provides a graphical warning (as in the case for the direct pathway from 1 → 2a) when important reaction intermediates have been missed. Although the NEB method is computationally intensive, it is also highly parallel, which is important in light of the fact that pathways connecting more than just the lowest energy adjacent minima should be explored. Insight generated from the application of the NEB method improved the computational modeling of the reaction pathway and holds significant promise for its ability to have a high impact in the field of rational catalyst design.
■
■
REFERENCES
(1) (a) Beller, M.; Bolm, C., Eds. Transition Metals for Organic Synthesis: Building Blocks and Fine Chemicals, Second Revised and Enlarged ed.; Wiley-VCH Verlag GmbH & Co. KGaA, 2004; Vol. 1; p 652. (b) Beller, M.; Bolm, C., Eds. Transition Metals for Organic Synthesis: Building Blocks and Fine Chemicals, Second ed.; Wiley-VCH Verlag GmbH & Co. KGaA, 2004; Vol. 2; p 662. (c) Ojima, I. Comprehensive Organometallic Chemistry III, Vol. 10: Applications II: Transition Metal Compounds in Organic Synthesis 1; Elsevier Ltd, 2007; p 891. (d) Schlosser, M.; Ed. Organometallics in Synthesis: A Manual, 2nd ed.; Wiley, 2002; p 1243. (e) Schlosser, M., Ed. Organometallics in Synthesis, Third Manual; John Wiley & Sons, Inc, 2013; p 1013. (f) Lipshutz, B. H., Ed. Organometallics in Synthesis, Fourth Manual; John Wiley & Sons, Inc: 2013; p 555. (2) Sperger, T.; Sanhueza, I. A.; Kalvet, I.; Schoenebeck, F. Chem. Rev. 2015, 115, 9532−9586. (3) (a) Hyatt, I. F. D.; Anderson, H. K.; Morehead, A. T.; Sargent, A. L. Organometallics 2008, 27, 135−147. (b) Delley, B. J. Chem. Phys. 1990, 92, 508−517. (c) Delley, B. J. Chem. Phys. 2000, 113, 7756−64. (d) Becke, A. D. J. Chem. Phys. 1988, 88, 1053−62. (e) Tsuneda, T.; Hirao, K. Chem. Phys. Lett. 1997, 268, 510−20. (f) Delley, B. Int. J. Quantum Chem. 1998, 69, 423−433. (4) (a) Fairlie, D. P.; Bosnich, B. Organometallics 1988, 7, 946−54. (b) Fairlie, D. P.; Bosnich, B. Organometallics 1988, 7, 936−45. (5) Alfonso, D. R.; Jordan, K. D. J. Comput. Chem. 2003, 24, 990−6. (6) (a) Jones, W. D. Comprehensive Organometallic Chemistry III, Vol. 1: 5. Advances in Carbon-Hydrogen Activation; Mingos, D. M. P.; Crabtree, R. H., Eds.; Elsevier Ltd, 2007; pp 699−723. (b) Hall, C.; Perutz, R. N. Chem. Rev. 1996, 96, 3125−46. (c) Ritleng, V.; Sirlin, C.; Pfeffer, M. Chem. Rev. 2002, 102, 1731−69. (7) (a) Sperger, T.; Sanhueza, I. A.; Kalvet, I.; Schoenebeck, F. Chem. Rev. 2015, 115, 9532−86. (b) Balcells, D.; Eisenstein, O. Comprehensive Inorganic Chemistry II; Poeppelmeier, J. R., Ed.; Elsevier: Amsterdam, 2013; Vol. 9, pp 695−726. (c) Niu, S.; Hall, M. B. Chem. Rev. 2000, 100, 353−405. (8) (a) Halgren, T. A.; Lipscomb, W. N. Chem. Phys. Lett. 1977, 49, 225−32. (b) Govind, N.; Peterson, M.; Fitzgerald, G.; King-Smith, D.; Andzelm, J. Comput. Mater. Sci. 2003, 28, 250−258. (9) (a) Jun, C.; Park, J. Top. Curr. Chem. 2013, 346, 59−84. (b) Jun, C.; Jo, E.; Park, J. Eur. J. Org. Chem. 2007, 2007, 1869−1881. (c) Lenges, C. P.; Brookhart, M. J. Am. Chem. Soc. 1999, 121, 6616− 23. (d) Moxham, G. L.; Randell-Sly, H. E.; Brayshaw, S. K.; Woodward, R. L.; Weller, A. S.; Willis, M. C. Angew. Chem., Int. Ed. 2006, 45, 7618−22. (e) Moxham, G.; Randell-Sly, H.; Brayshaw, S.; Weller, A.; Willis, M. Chem. - Eur. J. 2008, 14, 8383−97. (f) Osborne, J. D.; Willis, M. C. Chem. Commun. (Cambridge, U. K.) 2008, 5025− 27. (g) Park, Y. J.; Park, J.; Jun, C. Acc. Chem. Res. 2008, 41, 222−34. (h) Pike, S. D.; Pawley, R. J.; Chaplin, A. B.; Thompson, A. L.; Hooper, J. A.; Willis, M. C.; Weller, A. S. Eur. J. Inorg. Chem. 2011, 2011, 5558−65. (i) Prades, A.; Fernandez, M.; Pike, S. D.; Willis, M. C.; Weller, A. S. Angew. Chem., Int. Ed. 2015, 54, 8520−4. (j) Roy, A. H.; Lenges, C. P.; Brookhart, M. J. Am. Chem. Soc. 2007, 129, 2082− 2093. (k) Tanaka, K. Yuki Gosei Kagaku Kyokaishi 2012, 70, 1134−44. (l) Tanaka, K.; Fu, G. C. Org. Lett. 2002, 4, 933−935. (m) Chaplin, A. B.; Hooper, J. F.; Weller, A. S.; Willis, M. C. J. Am. Chem. Soc. 2012, 134, 4885−97. (10) The corresponding Gibbs free energies are listed in Figures S1 and S2 of the Supporting Information.
ASSOCIATED CONTENT
S Supporting Information *
The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.organomet.6b00236. The PES corresponding to Scheme 2 but with the Gibbs free energies for the oxidative addition pathways through 2a and 2b (PDF) Cartesian coordinates and total energies (in hartrees) of all structures presented herein using the standard xyz format (XYZ)
■
Article
AUTHOR INFORMATION
Corresponding Authors
*E-mail:
[email protected]. *E-mail:
[email protected]. Notes
The authors declare no competing financial interest. E
DOI: 10.1021/acs.organomet.6b00236 Organometallics XXXX, XXX, XXX−XXX
