Article pubs.acs.org/JPCA
Tunneling Above the Crossover Temperature Sonia Á lvarez-Barcia,† Jesús R. Flores,† and Johannes Kas̈ tner*,‡ †
Facultade de Química, Universidade de Vigo, E-36310 Vigo (Pontevedra), Spain Computational Biochemistry Group, Institute of Theoretical Chemistry, University of Stuttgart, Pfaffenwaldring 55, 70569 Stuttgart, Germany
‡
ABSTRACT: Quantum mechanical tunneling of atoms plays a significant role in many chemical reactions. The crossover temperature between classical and quantum movement is a convenient preliminary indication of the importance of tunneling for a particular reaction. Here we show, using instanton theory, that quantum tunneling is possible significantly above this crossover temperature for specific forms of the potential energy surface. We demonstrate the effect on an analytic potential as well as a chemical system. While protons move asynchronously along a Grotthuss chain in the classical high-temperature range, the onset of tunneling results in a synchronization of their movement.
■
INTRODUCTION Quantum tunneling, the transmission through energy barriers higher than the energy of the involved particles, plays a role in many chemical reactions.1,2 Tunneling competes with the classical thermally activated barrier-crossing. Since tunneling is temperature independent, it is most important at low temperature. However, many cases of tunneling at room temperature in chemical and biochemical systems have been described3,4 as well as processes at cryogenic temperature.5,6 The theoretical prediction of tunneling rates is much more elaborate than the prediction of classical rates. The most general and rigorous approach is, of course, to perform quantum dynamics, i.e., to integrate the time-dependent Schrödinger equation starting from the reactant state ensemble. Due to the computational demand, this approach is usually restricted to molecular systems with only very few atoms. Alternatives are the Feynman-path-based methods centroid molecular dynamics7,8 or ring-polymer molecular dynamics.9 The semiclassical approximation often provides sufficient accuracy to describe the tunneling motion of atoms. It focuses on one main tunneling path. If that path is chosen close to the classical intrinsic reaction coordinate, the small-curvature tunneling correction (SCT)10 results. If it is chosen as the straight-line path, minimizing the barrier width, the resulting method is called large-curvature tunneling correction (LCT).11,12 The tunneling path is optimized for a maximum in the tunneling probability by the instanton method,13−17 which is used in the present work. It is based on Feynman path integral theory and employs the semiclassical approximation. Before any of these methods are employed, however, it is advisable to estimate the importance of tunneling for a specific reaction. The crossover temperature Tc may be employed for this task:18 © 2013 American Chemical Society
Tc =
ℏΩ 2πkB
(1)
with Ω being the barrier frequency (the absolute value of the imaginary frequency corresponding to the transition mode) and kB corresponding to Boltzmann’s constant. Calculating Tc, thus, requires only quantities which are known whenever a transition state is characterized. Tc generally marks the temperature where tunneling and thermally activated barrier crossing are roughly equally important. As an example, Tc is roughly room temperature for reactions with Ω = 1300 cm−1. Using this analysis, a significant tunneling contribution is at room temperature expected only for transitions with a higher barrier frequency. The method of choice for the present work to calculate tunneling rates is instanton theory.13−17 In a nutshell, statistical Feynman path integral theory is used together with transition state theory to describe quantum ensembles of the reactant state and the transition state. The instanton is the most probable tunneling path. For more details, the reader is referred to the original literature13−17 or some recommended derivative work.2,19−21 Instanton theory is applicable whenever the temperature is low enough for the instanton to spread out. At high enough temperature, the instanton collapses to a point, which renders the theory inapplicable. For many barrier shapes this collapse happens at the crossover temperature Tc. Here, we will show examples in which this happens, but more importantly, also examples where delocalized instanton solutions are possible above Tc. The paper is organized as follows. In the next section, analytic examples are discussed which show under which Received: November 14, 2013 Published: December 12, 2013 78
dx.doi.org/10.1021/jp411189m | J. Phys. Chem. A 2014, 118, 78−82
The Journal of Physical Chemistry A
Article
conditions stable instanton solutions are possible above Tc. Then a chemical system for which this is indeed the case is introduced, the hydration of aluminum by water via a Grotthuss chain.22−27
■
METHODS Theory. The instanton is a first-order saddle point of the Euclidean action13−17 SE =
ℏ kBT
P
⎛ P(k T )2 B
∑ ⎜⎜ k=1
⎝
2
2ℏ
|yk + 1 − yk |2 +
V (yk ) ⎞ ⎟⎟ P ⎠
(2)
with kB being Boltzmann’s constant and T being the absolute temperature, ℏ is Planck’s constant, P is the number of images to discretize the Feynman path (instanton path), yk are the mass-weighted coordinates of image k of the Feynman path, and V(y) is the potential energy. Mass weighted coordinates are obtained by yi = (mi)1/2xi from the Cartesian coordinates xi and the mass mi of the atom to which component i belongs. The lower SE the higher is the tunneling probability. The first term in the parentheses of eq 2, (P(kBT)2/2ℏ2)|yk+1 − yk|2, causes the instanton path to shrink. Its effect is stronger the higher the temperature is. The second term, the potential energy, causes the path to expand. The instanton path stretches over the barrier, its ends are both at the energy Eb, the tunneling energy. The longer the path, the more images have a low potential energy. The second term is temperature independent. At some temperature, the first term dominates over the second and a short path collapses to a point, a saddle point on the potential energy surface. If the potential drops less steeply than the harmonic approximation at this saddle point, the path collapses at the crossover temperature Tc. Since this is the case for most chemical reactions, Tc is often a good estimate for the temperature at which tunneling starts to play a dominant role. Analytic Examples. Here, analytic examples with onedimensional reaction coordinates are discussed. While in real chemical systems, each atom has three degrees of freedom, which are all accounted for, these examples may be regarded as effective one-dimensional Hamiltonians. A one-dimensional example in which the instanton collapses to the saddle point at T = Tc is V(y) = −2y3 − 3y2. It is shown in Figure 1. All quantities are given in atomic units. A mass of the hydrogen atom was used. Toward the reactant state (left), V(y) is larger and less steep than the quadratic approximation. The dependence of the tunneling energy Eb on T is also shown in Figure 1, top right panel. For T > Tc the instanton is collapsed and Eb is formally equal to the energy of the barrier top (V = 0 in this case). At lower temperature, Eb drops and finally reaches the energy of the reactant state minium (V = −1 in this case). The reaction rate, shown as an Arrhenius plot depicting the logarithm of the rate plotted vs the inverse temperature, is shown on the bottom right of Figure 1. One can, however, also imagine cases in which the potential is rather flat at the top of a barrier, but decreases steeper than the quadratic approximation further away from the barrier. The analytic function V(y) = −2y5 − 2.5y4 − 0.5y2 may serve as an example. It is illustrated in Figure 2. The potential decreases steeper than the harmonic approximation from the transition state toward both reactant and product states. In this case one can find delocalized instanton solutions even for temperatures above Tc. The temperature dependence of the rate is depicted on the lower (right) panel of Figure 2. It looks qualitatively
Figure 1. Top left: Energy vs reaction coordinate y of the potential V(y) = −2y3 − 3y2. Top right: Tunneling energy Eb vs inverse temperature (Tc/T). Bottom right: Arrhenius plot, ln(rate) vs Tc/T. This is a typical behavior of an instanton.
Figure 2. Top left: Energy vs reaction coordinate y of the potential V(y) = −2y5 − 2.5y4 − 0.5y2. Top right: Tunneling energy Eb vs inverse temperature (Tc/T). Bottom right: Arrhenius plot, ln(rate) vs Tc/T. Here, tunneling is observed above Tc.
similar to the usual case of Figure 1, with the only difference that tunneling takes place also for T > Tc. Only at about T = 1.82Tc does the instanton solution collapse. To discuss the effects in this potential it makes sense to start at the low-temperature limit. As described above, the instanton spreads from the reactant minimum to the product side. At low temperature, Eb is very close to the energy of the reactant. With higher temperature, Eb increases. For potentials of the shape shown in Figure 2, though, two possible stationary instanton paths for the same temperature, but very different Eb can be found for Tc < T < 1.82Tc. The solution for the higher Eb is a second-order saddle point, though, and either collapses to a point (at the classical saddle point) or spreads to the other solution, the one with the lower Eb. The eigenvalue spectrum of the unstable solution is wrong, another negative eigenvalue of the Hessian is found, which corresponds to a spreading of the instanton path. The solution with the lower Eb is a valid instanton with a correct eigenvalue spectrum. It results in tunneling even above the crossover temperature. The present case is not a limitation of instanton theory. For such potentials tunneling will, in fact, occur at T > Tc. It is more a hint that caution is needed when interpreting Tc. 79
dx.doi.org/10.1021/jp411189m | J. Phys. Chem. A 2014, 118, 78−82
The Journal of Physical Chemistry A
■
Article
COMPUTATIONAL DETAILS All geometry optimizations and instanton optimizations were performed with DL-FIND20,28 in ChemShell.29−31 The calculations were carried out at the M06/6-311+G** level.32−34 For the electronic structure computations we have employed the G09 program.35 SCF cycles were stopped when the convergence, as defined in G09,35 reached 10−8. We have also used a pruned (99 590) grid, having 99 radial shells and 590 angular points per shell.
With lowering of the temperature, the part of the instanton path with the highest energy deviates more and more from the IRC. This leads to a higher, but thinner barrier for instantons at low temperature, see Figure 4. A thin barrier has a higher tunneling probability than a wide one. The potential energy barrier for a proton transfer from O3 via O2 and O1 to aluminum is 12.6 kcal mol−1. The highest points on the temperature-dependent tunneling paths are 13.8 to 18.2 kcal mol−1 above the energy of the reactant, depending on the temperature. Thus, even such an increase in the energy is counterbalanced by a reduction in the barrier thickness. The geometric changes between the instanton path and the classical transition state can be seen in Figure 3 on the left panel. It shows the instanton at the highest temperature that was found to be stable in delocalized form compared to the classical transition state. In the instanton tunneling path, all protons are transferred in a synchronized manner. The path is shown starting in blue and ending in red. In the classical reaction, however, first the proton between O3 and O2 is transferred, while the one between O1 and Al is still bound to O1. Only later during the mechanism is the proton finally moved to Al. Thus, tunneling turns an asynchronous protonshift reaction along a Grotthuss chain into a synchronous mechanism. This occurs at the expense of a significant amount of energy (1.2 to 5.6 kcal mol−1), but leads to a reduction of the barrier widths. The tunneling reaction takes a synchronous path in the relay chain, even at the highest temperature for which instantons could be localized. All hydrogen atoms are transferred synchronously. The difference between the two paths can also be depicted geometrically, see Figure 5. The tunneling
■
RESULTS Tunneling above Tc, described by delocalized instanton solutions, does not only occur in analytic potentials as shown in the Analytic Examples section, but can also be found in real molecular systems. An example is shown here. It has been found36 that the hydrogen transfer reaction in microhydrated aluminum atoms can proceed very fast by a relay (Grotthuss-like) reaction mechanism, as depicted in Figure 3.
Figure 3. Al + 3H2O: instanton at 230.5 K compared to the classical TS (white, positions of the H atoms are indicated by black dashes). Note that in the classical transition state, the hydrogen atom between O2 and O3 is almost transferred to O2, while the one between O1 and Al still remains close to O1.
Hydrogen transfer from water to aluminum was modeled by one aluminum atom and three water molecules, see Figure 3. In this multidimensional system, a crossover temperature of Tc = 149.3 K was obtained, while delocalized instantons were found up to 230 K. This is caused by two effects, both can be seen in Figure 4. First, the potential along the intrinsic reaction coordinate (IRC) decreases steeper than the harmonic approximation. This is the same as observed in the second analytic example above. However, in a multidimensional system like this one, the instanton path can also qualitatively deviate from the IRC. This is the case here. Figure 5. Corner cutting in the instanton path for Al + 3H2O: the instantons do not proceed through the classical transition state (circle) but shorten the path toward the product state (upper right) in a more synchronous mechanism. Stationary points are indicated by circles.
paths for two different temperatures are compared to the IRC projected to the difference of distances d(O1−H) − d(Al−H) and d(O3−H) − d(O2−H). These coordinates describe the hydrogen transfer reactions from O1 to Al and from O3 to O2, respectively. It can clearly be seen that for both temperatures depicted, tunneling proceeds via a different, more synchronized, route than the classical asynchronous reaction. The resulting rate constants are shown in Figure 6. Tunneling sets in at about 230 K, slightly below room temperature. With lowering temperature, the rate constant drops by about 4 orders of magnitude to 60 s−1 or a half-life of
Figure 4. Energy along the tunneling paths at different temperatures. Higher barriers are crossed at lower temperature. The IRC and the harmonic approximation around the transition state are indicated as well. 80
dx.doi.org/10.1021/jp411189m | J. Phys. Chem. A 2014, 118, 78−82
The Journal of Physical Chemistry A
Article
Excellence in Simulation Technology (EXC 310/1) at the University of Stuttgart.
■
(1) Ley, D.; Gerbig, D.; Schreiner, P. R. Tunnelling control of chemical reactions − the organic chemist’s perspective. Org. Biomol. Chem. 2012, 10, 3781. (2) Kästner, J. Theory and Simulation of Atom Tunneling in Chemical Reactions. Comput. Mol. Sci. 2013, DOI: 10.1002/ wcms.1165. (3) Sutcliffe, M. J.; Scrutton, N. S. A new conceptual framework for enzyme catalysis − hydrogen tunneling coupled to enzyme dynamics in flavoprotein and quinoprotein enzymes. Eur. J. Biochem. 2002, 269, 3096. (4) Masgrau, L.; Roujeinikova, A.; Johannissen, L. O.; Hothi, P.; Basran, J.; Ranaghan, K. E.; Mulholland, A. J.; Sutcliffe, M. J.; Scrutton, N. S.; Leys, D. Atomic Description of an Enzyme Reaction Dominated by Proton Tunneling. Science 2006, 312, 237. (5) Schreiner, P. R.; Reisenauer, H. P.; Pickard, F. C.; Simmonett, A. C.; Allen, W. D.; Mátyus, E.; Császár, A. G. Capture of hydroxymethylene and its fast disappearance through tunnelling. Nature 2008, 453, 906. (6) Schreiner, P. R.; Reisenauer, H. P.; Ley, D.; Gerbig, D.; Wu, C.H.; Allen, W. D. Methylhydroxycarbene: Tunneling Control of a Chemical Reaction. Science 2011, 332, 1300. (7) Cao, J.; Voth, G. A. The formulation of quantum statistical mechanics based on the Feynman path centroid density. III. Phase space formalism and analysis of centroid molecular dynamics. J. Chem. Phys. 1994, 101, 6157. (8) Voth, G. A. Path-Integral Centroid Methods in Quantum Statistical Mechanics and Dynamics. Adv. Chem. Phys. 1996, 93, 135. (9) Craig, I. R.; Manolopoulos, D. E. Chemical reaction rates from ring polymer molecular dynamics. J. Chem. Phys. 2005, 122, 084106. (10) Skodje, R. T.; Truhlar, D. G.; Garrett, B. C. A general smallcurvature approximation for transition-state-theory transmission coefficients. J. Phys. Chem. 1981, 85, 3019−3023. (11) Garrett, B. C.; Truhlar, D. G.; Wagner, A. F.; Dunning, T. H., Jr. Variational transition state theory and tunneling for a heavy−light− heavy reaction using an ab initio potential energy surface. J. Chem. Phys. 1983, 78, 4400. (12) Garrett, B. C.; Abusalbi, N.; Kouri, D. J.; Truhlar, D. G. Test of variational transition state theory and the least-action approximation for multidimensional tunneling probabilities against accurate quantal rate constants for a collinear reaction involving tunneling into an excited state. J. Chem. Phys. 1985, 83, 2252. (13) Langer, J. S. Theory of the condensation point. Ann. Phys. (N.Y.) 1967, 41, 108. (14) Miller, W. H. Semiclassical limit of quantum mechanical transition state theory for nonseparable systems. J. Chem. Phys. 1975, 62, 1899. (15) Coleman, S. Fate of the false vacuum: Semiclassical theory. Phys. Rev. D 1977, 15, 2929. (16) Callan, C. G., Jr.; Coleman, S. Fate of the false vacuum. II. First quantum corrections. Phys. Rev. D 1977, 16, 1762. (17) Gildener, E.; Patrascioiu, A. Pseudoparticle contributions to the energy spectrum of a one-dimensional system. Phys. Rev. D 1977, 16, 423. (18) Gillan, M. J. Quantum-classical crossover of the transition rate in the damped double well. J. Phys. C 1987, 20, 3621. (19) Richardson, J. O.; Althorpe, S. C. Ring-polymer molecular dynamics rate-theory in the deep-tunneling regime: Connection with semiclassical instanton theory. J. Chem. Phys. 2009, 131, 214106. (20) Rommel, J. B.; Goumans, T. P. M.; Kästner, J. Locating instantons in many degrees of freedom. J. Chem. Theory Comput. 2011, 7, 690. (21) Rommel, J. B.; Kästner, J. Adaptive Integration Grids in Instanton Theory Improve the Numerical Accuracy at Low Temperature. J. Chem. Phys. 2011, 134, 184107.
Figure 6. Arrhenius plot of the rate constant vs inverse temperature.
11.5 ms at a temperature of 50 K. The classical rate constant (which includes zero point vibrational energy, but no tunneling contributions) is many orders of magnitude smaller than the tunneling rate constant. The remaining slight temperature dependence of the tunneling rate constant shows that, while tunneling leads to synchronous concerted proton transfers, it is still aided by the finite thermal motion, a process known as temperature-assisted tunneling. The physical background is that in the reactant state, all molecules are further apart from each other than energetically favorable for the proton transfers. The contraction of the complex requires hardly any energy, but significant movements of atoms. This process is more efficient in classical movement than in tunneling of atoms down to very low temperature. The effect shows in long flat tails in the energy curves of the instantons at low temperature to the left in Figure 4 and in the still noticeable decrease of the tunneling rate at very low temperature.
■
CONCLUSION We have shown that quantum tunneling of atoms synchronizes the proton movements along Grotthuss chains. While the chosen example is a special case of a strongly exothermic reaction with an aluminum atom at the end of the transport chain, a similar mechanism may be relevant in liquid water. We have used instanton theory to calculate tunneling rate constants. Extended instanton solutions and corresponding tunneling rates were obtained in analytic cases and in a molecular example even above the crossover temperature Tc. While Tc is still a helpful concept to judge without expensive calculations if tunneling might play a role in a particular system, the present study shows that such analysis has to be interpreted with care.
■
REFERENCES
AUTHOR INFORMATION
Corresponding Author
*E-mail: [email protected]. Phone: +49 711 685 64473. Fax: +49 711 685 64442. Notes
The authors declare no competing financial interest.
■
ACKNOWLEDGMENTS J.R.F. and S.A.B. acknowledge the financial support of the University of Vigo (Contrato-Programa, 2012). S.A.B. acknowledges a F.P.U. grant from the Spanish Ministry of Education. J.K. would like to thank the German Research Foundation (DFG) for financial support of the project within the Cluster of 81
dx.doi.org/10.1021/jp411189m | J. Phys. Chem. A 2014, 118, 78−82
The Journal of Physical Chemistry A
Article
(22) de Grotthuss, C. J. T. Sur la décomposition de l’eau et des corps qu’elle tient en dissolution à l’aide de l’électricité galvanique. Ann. Chim. (Paris) 1806, 58, 54. (23) Hückel, E. Theory of the mobilities of the hydrogen and hydroxyl ions in aqueous solution. Z. Electrochem. 1928, 34, 546. (24) Bernal, J. D.; Fowler, R. H. A Theory of Water and Ionic Solution, with Particular Reference to Hydrogen and Hydroxyl Ions. J. Chem. Phys. 1933, 1, 515. (25) Agmon, N. The Grotthuss mechanism. Chem. Phys. Lett. 1995, 244, 456. (26) Cukierman, S. Et tu, Grotthuss! and other unfinished stories. Biochim. Biophys. Acta 2006, 1757, 876. (27) Knight, C.; Voth, G. A. The Curious Case of the Hydrated Proton. Acc. Chem. Res. 2012, 45, 101. (28) Kästner, J.; Carr, J. M.; Keal, T. W.; Thiel, W.; Wander, A.; Sherwood, P. DL-FIND: an Open-Source Geometry Optimizer for Atomistic Simulations. J. Phys. Chem. A 2009, 113, 11856. (29) Sherwood, P.; et al. QUASI: A general purpose implementation of the QM/MM approach and its application to problems in catalysis. J. Mol. Struct. (THEOCHEM) 2003, 632, 1. (30) Metz, S.; Kästner, J.; Sokol, A. A.; Keal, T. W.; Sherwood, P. ChemShell−a modular software package for QM/MM simulations. Comput. Mol. Sci. 2013, DOI: 10.1002/wcms.1163. (31) ChemShell, a Computational Chemistry Shell, see www. chemshell.org. (32) Zhao, Y.; Truhlar, D. G. The M06 suite of density functionals for main group thermochemistry, thermochemical kinetics, noncovalent interactions, excited states, and transition elements: two new functionals and systematic testing of four M06-class functionals and 12 other functionals. Theor. Chem. Acc. 2008, 120, 215. (33) Clark, T.; Chandrasekhar, J.; Spitznagel, G. W.; Schleyer, P. v. R. Efficient diffuse function-augmented basis sets for anion calculations. III. The 3-21+G basis set for first-row elements, Li-F. J. Comput. Chem. 1983, 4, 294. (34) Krishnan, R.; Binkley, J. S.; Seeger, R.; Pople, J. A. Selfconsistent molecular orbital methods. XX. A basis set for correlated wave functions. J. Chem. Phys. 1980, 72, 650. (35) Frisch, M. J. et al. Gaussian 09, Revision C.01. (36) Á lvarez-Barcia, S.; Flores, J. R. The interaction of Al atoms with water molecules: A theoretical study. J. Chem. Phys. 2009, 131, 174307.
82
dx.doi.org/10.1021/jp411189m | J. Phys. Chem. A 2014, 118, 78−82
Tunneling Above the Crossover Temperature Sonia Á lvarez-Barcia,† Jesús R. Flores,† and Johannes Kas̈ tner*,‡ †
Facultade de Química, Universidade de Vigo, E-36310 Vigo (Pontevedra), Spain Computational Biochemistry Group, Institute of Theoretical Chemistry, University of Stuttgart, Pfaffenwaldring 55, 70569 Stuttgart, Germany
‡
ABSTRACT: Quantum mechanical tunneling of atoms plays a significant role in many chemical reactions. The crossover temperature between classical and quantum movement is a convenient preliminary indication of the importance of tunneling for a particular reaction. Here we show, using instanton theory, that quantum tunneling is possible significantly above this crossover temperature for specific forms of the potential energy surface. We demonstrate the effect on an analytic potential as well as a chemical system. While protons move asynchronously along a Grotthuss chain in the classical high-temperature range, the onset of tunneling results in a synchronization of their movement.
■
INTRODUCTION Quantum tunneling, the transmission through energy barriers higher than the energy of the involved particles, plays a role in many chemical reactions.1,2 Tunneling competes with the classical thermally activated barrier-crossing. Since tunneling is temperature independent, it is most important at low temperature. However, many cases of tunneling at room temperature in chemical and biochemical systems have been described3,4 as well as processes at cryogenic temperature.5,6 The theoretical prediction of tunneling rates is much more elaborate than the prediction of classical rates. The most general and rigorous approach is, of course, to perform quantum dynamics, i.e., to integrate the time-dependent Schrödinger equation starting from the reactant state ensemble. Due to the computational demand, this approach is usually restricted to molecular systems with only very few atoms. Alternatives are the Feynman-path-based methods centroid molecular dynamics7,8 or ring-polymer molecular dynamics.9 The semiclassical approximation often provides sufficient accuracy to describe the tunneling motion of atoms. It focuses on one main tunneling path. If that path is chosen close to the classical intrinsic reaction coordinate, the small-curvature tunneling correction (SCT)10 results. If it is chosen as the straight-line path, minimizing the barrier width, the resulting method is called large-curvature tunneling correction (LCT).11,12 The tunneling path is optimized for a maximum in the tunneling probability by the instanton method,13−17 which is used in the present work. It is based on Feynman path integral theory and employs the semiclassical approximation. Before any of these methods are employed, however, it is advisable to estimate the importance of tunneling for a specific reaction. The crossover temperature Tc may be employed for this task:18 © 2013 American Chemical Society
Tc =
ℏΩ 2πkB
(1)
with Ω being the barrier frequency (the absolute value of the imaginary frequency corresponding to the transition mode) and kB corresponding to Boltzmann’s constant. Calculating Tc, thus, requires only quantities which are known whenever a transition state is characterized. Tc generally marks the temperature where tunneling and thermally activated barrier crossing are roughly equally important. As an example, Tc is roughly room temperature for reactions with Ω = 1300 cm−1. Using this analysis, a significant tunneling contribution is at room temperature expected only for transitions with a higher barrier frequency. The method of choice for the present work to calculate tunneling rates is instanton theory.13−17 In a nutshell, statistical Feynman path integral theory is used together with transition state theory to describe quantum ensembles of the reactant state and the transition state. The instanton is the most probable tunneling path. For more details, the reader is referred to the original literature13−17 or some recommended derivative work.2,19−21 Instanton theory is applicable whenever the temperature is low enough for the instanton to spread out. At high enough temperature, the instanton collapses to a point, which renders the theory inapplicable. For many barrier shapes this collapse happens at the crossover temperature Tc. Here, we will show examples in which this happens, but more importantly, also examples where delocalized instanton solutions are possible above Tc. The paper is organized as follows. In the next section, analytic examples are discussed which show under which Received: November 14, 2013 Published: December 12, 2013 78
dx.doi.org/10.1021/jp411189m | J. Phys. Chem. A 2014, 118, 78−82
The Journal of Physical Chemistry A
Article
conditions stable instanton solutions are possible above Tc. Then a chemical system for which this is indeed the case is introduced, the hydration of aluminum by water via a Grotthuss chain.22−27
■
METHODS Theory. The instanton is a first-order saddle point of the Euclidean action13−17 SE =
ℏ kBT
P
⎛ P(k T )2 B
∑ ⎜⎜ k=1
⎝
2
2ℏ
|yk + 1 − yk |2 +
V (yk ) ⎞ ⎟⎟ P ⎠
(2)
with kB being Boltzmann’s constant and T being the absolute temperature, ℏ is Planck’s constant, P is the number of images to discretize the Feynman path (instanton path), yk are the mass-weighted coordinates of image k of the Feynman path, and V(y) is the potential energy. Mass weighted coordinates are obtained by yi = (mi)1/2xi from the Cartesian coordinates xi and the mass mi of the atom to which component i belongs. The lower SE the higher is the tunneling probability. The first term in the parentheses of eq 2, (P(kBT)2/2ℏ2)|yk+1 − yk|2, causes the instanton path to shrink. Its effect is stronger the higher the temperature is. The second term, the potential energy, causes the path to expand. The instanton path stretches over the barrier, its ends are both at the energy Eb, the tunneling energy. The longer the path, the more images have a low potential energy. The second term is temperature independent. At some temperature, the first term dominates over the second and a short path collapses to a point, a saddle point on the potential energy surface. If the potential drops less steeply than the harmonic approximation at this saddle point, the path collapses at the crossover temperature Tc. Since this is the case for most chemical reactions, Tc is often a good estimate for the temperature at which tunneling starts to play a dominant role. Analytic Examples. Here, analytic examples with onedimensional reaction coordinates are discussed. While in real chemical systems, each atom has three degrees of freedom, which are all accounted for, these examples may be regarded as effective one-dimensional Hamiltonians. A one-dimensional example in which the instanton collapses to the saddle point at T = Tc is V(y) = −2y3 − 3y2. It is shown in Figure 1. All quantities are given in atomic units. A mass of the hydrogen atom was used. Toward the reactant state (left), V(y) is larger and less steep than the quadratic approximation. The dependence of the tunneling energy Eb on T is also shown in Figure 1, top right panel. For T > Tc the instanton is collapsed and Eb is formally equal to the energy of the barrier top (V = 0 in this case). At lower temperature, Eb drops and finally reaches the energy of the reactant state minium (V = −1 in this case). The reaction rate, shown as an Arrhenius plot depicting the logarithm of the rate plotted vs the inverse temperature, is shown on the bottom right of Figure 1. One can, however, also imagine cases in which the potential is rather flat at the top of a barrier, but decreases steeper than the quadratic approximation further away from the barrier. The analytic function V(y) = −2y5 − 2.5y4 − 0.5y2 may serve as an example. It is illustrated in Figure 2. The potential decreases steeper than the harmonic approximation from the transition state toward both reactant and product states. In this case one can find delocalized instanton solutions even for temperatures above Tc. The temperature dependence of the rate is depicted on the lower (right) panel of Figure 2. It looks qualitatively
Figure 1. Top left: Energy vs reaction coordinate y of the potential V(y) = −2y3 − 3y2. Top right: Tunneling energy Eb vs inverse temperature (Tc/T). Bottom right: Arrhenius plot, ln(rate) vs Tc/T. This is a typical behavior of an instanton.
Figure 2. Top left: Energy vs reaction coordinate y of the potential V(y) = −2y5 − 2.5y4 − 0.5y2. Top right: Tunneling energy Eb vs inverse temperature (Tc/T). Bottom right: Arrhenius plot, ln(rate) vs Tc/T. Here, tunneling is observed above Tc.
similar to the usual case of Figure 1, with the only difference that tunneling takes place also for T > Tc. Only at about T = 1.82Tc does the instanton solution collapse. To discuss the effects in this potential it makes sense to start at the low-temperature limit. As described above, the instanton spreads from the reactant minimum to the product side. At low temperature, Eb is very close to the energy of the reactant. With higher temperature, Eb increases. For potentials of the shape shown in Figure 2, though, two possible stationary instanton paths for the same temperature, but very different Eb can be found for Tc < T < 1.82Tc. The solution for the higher Eb is a second-order saddle point, though, and either collapses to a point (at the classical saddle point) or spreads to the other solution, the one with the lower Eb. The eigenvalue spectrum of the unstable solution is wrong, another negative eigenvalue of the Hessian is found, which corresponds to a spreading of the instanton path. The solution with the lower Eb is a valid instanton with a correct eigenvalue spectrum. It results in tunneling even above the crossover temperature. The present case is not a limitation of instanton theory. For such potentials tunneling will, in fact, occur at T > Tc. It is more a hint that caution is needed when interpreting Tc. 79
dx.doi.org/10.1021/jp411189m | J. Phys. Chem. A 2014, 118, 78−82
The Journal of Physical Chemistry A
■
Article
COMPUTATIONAL DETAILS All geometry optimizations and instanton optimizations were performed with DL-FIND20,28 in ChemShell.29−31 The calculations were carried out at the M06/6-311+G** level.32−34 For the electronic structure computations we have employed the G09 program.35 SCF cycles were stopped when the convergence, as defined in G09,35 reached 10−8. We have also used a pruned (99 590) grid, having 99 radial shells and 590 angular points per shell.
With lowering of the temperature, the part of the instanton path with the highest energy deviates more and more from the IRC. This leads to a higher, but thinner barrier for instantons at low temperature, see Figure 4. A thin barrier has a higher tunneling probability than a wide one. The potential energy barrier for a proton transfer from O3 via O2 and O1 to aluminum is 12.6 kcal mol−1. The highest points on the temperature-dependent tunneling paths are 13.8 to 18.2 kcal mol−1 above the energy of the reactant, depending on the temperature. Thus, even such an increase in the energy is counterbalanced by a reduction in the barrier thickness. The geometric changes between the instanton path and the classical transition state can be seen in Figure 3 on the left panel. It shows the instanton at the highest temperature that was found to be stable in delocalized form compared to the classical transition state. In the instanton tunneling path, all protons are transferred in a synchronized manner. The path is shown starting in blue and ending in red. In the classical reaction, however, first the proton between O3 and O2 is transferred, while the one between O1 and Al is still bound to O1. Only later during the mechanism is the proton finally moved to Al. Thus, tunneling turns an asynchronous protonshift reaction along a Grotthuss chain into a synchronous mechanism. This occurs at the expense of a significant amount of energy (1.2 to 5.6 kcal mol−1), but leads to a reduction of the barrier widths. The tunneling reaction takes a synchronous path in the relay chain, even at the highest temperature for which instantons could be localized. All hydrogen atoms are transferred synchronously. The difference between the two paths can also be depicted geometrically, see Figure 5. The tunneling
■
RESULTS Tunneling above Tc, described by delocalized instanton solutions, does not only occur in analytic potentials as shown in the Analytic Examples section, but can also be found in real molecular systems. An example is shown here. It has been found36 that the hydrogen transfer reaction in microhydrated aluminum atoms can proceed very fast by a relay (Grotthuss-like) reaction mechanism, as depicted in Figure 3.
Figure 3. Al + 3H2O: instanton at 230.5 K compared to the classical TS (white, positions of the H atoms are indicated by black dashes). Note that in the classical transition state, the hydrogen atom between O2 and O3 is almost transferred to O2, while the one between O1 and Al still remains close to O1.
Hydrogen transfer from water to aluminum was modeled by one aluminum atom and three water molecules, see Figure 3. In this multidimensional system, a crossover temperature of Tc = 149.3 K was obtained, while delocalized instantons were found up to 230 K. This is caused by two effects, both can be seen in Figure 4. First, the potential along the intrinsic reaction coordinate (IRC) decreases steeper than the harmonic approximation. This is the same as observed in the second analytic example above. However, in a multidimensional system like this one, the instanton path can also qualitatively deviate from the IRC. This is the case here. Figure 5. Corner cutting in the instanton path for Al + 3H2O: the instantons do not proceed through the classical transition state (circle) but shorten the path toward the product state (upper right) in a more synchronous mechanism. Stationary points are indicated by circles.
paths for two different temperatures are compared to the IRC projected to the difference of distances d(O1−H) − d(Al−H) and d(O3−H) − d(O2−H). These coordinates describe the hydrogen transfer reactions from O1 to Al and from O3 to O2, respectively. It can clearly be seen that for both temperatures depicted, tunneling proceeds via a different, more synchronized, route than the classical asynchronous reaction. The resulting rate constants are shown in Figure 6. Tunneling sets in at about 230 K, slightly below room temperature. With lowering temperature, the rate constant drops by about 4 orders of magnitude to 60 s−1 or a half-life of
Figure 4. Energy along the tunneling paths at different temperatures. Higher barriers are crossed at lower temperature. The IRC and the harmonic approximation around the transition state are indicated as well. 80
dx.doi.org/10.1021/jp411189m | J. Phys. Chem. A 2014, 118, 78−82
The Journal of Physical Chemistry A
Article
Excellence in Simulation Technology (EXC 310/1) at the University of Stuttgart.
■
(1) Ley, D.; Gerbig, D.; Schreiner, P. R. Tunnelling control of chemical reactions − the organic chemist’s perspective. Org. Biomol. Chem. 2012, 10, 3781. (2) Kästner, J. Theory and Simulation of Atom Tunneling in Chemical Reactions. Comput. Mol. Sci. 2013, DOI: 10.1002/ wcms.1165. (3) Sutcliffe, M. J.; Scrutton, N. S. A new conceptual framework for enzyme catalysis − hydrogen tunneling coupled to enzyme dynamics in flavoprotein and quinoprotein enzymes. Eur. J. Biochem. 2002, 269, 3096. (4) Masgrau, L.; Roujeinikova, A.; Johannissen, L. O.; Hothi, P.; Basran, J.; Ranaghan, K. E.; Mulholland, A. J.; Sutcliffe, M. J.; Scrutton, N. S.; Leys, D. Atomic Description of an Enzyme Reaction Dominated by Proton Tunneling. Science 2006, 312, 237. (5) Schreiner, P. R.; Reisenauer, H. P.; Pickard, F. C.; Simmonett, A. C.; Allen, W. D.; Mátyus, E.; Császár, A. G. Capture of hydroxymethylene and its fast disappearance through tunnelling. Nature 2008, 453, 906. (6) Schreiner, P. R.; Reisenauer, H. P.; Ley, D.; Gerbig, D.; Wu, C.H.; Allen, W. D. Methylhydroxycarbene: Tunneling Control of a Chemical Reaction. Science 2011, 332, 1300. (7) Cao, J.; Voth, G. A. The formulation of quantum statistical mechanics based on the Feynman path centroid density. III. Phase space formalism and analysis of centroid molecular dynamics. J. Chem. Phys. 1994, 101, 6157. (8) Voth, G. A. Path-Integral Centroid Methods in Quantum Statistical Mechanics and Dynamics. Adv. Chem. Phys. 1996, 93, 135. (9) Craig, I. R.; Manolopoulos, D. E. Chemical reaction rates from ring polymer molecular dynamics. J. Chem. Phys. 2005, 122, 084106. (10) Skodje, R. T.; Truhlar, D. G.; Garrett, B. C. A general smallcurvature approximation for transition-state-theory transmission coefficients. J. Phys. Chem. 1981, 85, 3019−3023. (11) Garrett, B. C.; Truhlar, D. G.; Wagner, A. F.; Dunning, T. H., Jr. Variational transition state theory and tunneling for a heavy−light− heavy reaction using an ab initio potential energy surface. J. Chem. Phys. 1983, 78, 4400. (12) Garrett, B. C.; Abusalbi, N.; Kouri, D. J.; Truhlar, D. G. Test of variational transition state theory and the least-action approximation for multidimensional tunneling probabilities against accurate quantal rate constants for a collinear reaction involving tunneling into an excited state. J. Chem. Phys. 1985, 83, 2252. (13) Langer, J. S. Theory of the condensation point. Ann. Phys. (N.Y.) 1967, 41, 108. (14) Miller, W. H. Semiclassical limit of quantum mechanical transition state theory for nonseparable systems. J. Chem. Phys. 1975, 62, 1899. (15) Coleman, S. Fate of the false vacuum: Semiclassical theory. Phys. Rev. D 1977, 15, 2929. (16) Callan, C. G., Jr.; Coleman, S. Fate of the false vacuum. II. First quantum corrections. Phys. Rev. D 1977, 16, 1762. (17) Gildener, E.; Patrascioiu, A. Pseudoparticle contributions to the energy spectrum of a one-dimensional system. Phys. Rev. D 1977, 16, 423. (18) Gillan, M. J. Quantum-classical crossover of the transition rate in the damped double well. J. Phys. C 1987, 20, 3621. (19) Richardson, J. O.; Althorpe, S. C. Ring-polymer molecular dynamics rate-theory in the deep-tunneling regime: Connection with semiclassical instanton theory. J. Chem. Phys. 2009, 131, 214106. (20) Rommel, J. B.; Goumans, T. P. M.; Kästner, J. Locating instantons in many degrees of freedom. J. Chem. Theory Comput. 2011, 7, 690. (21) Rommel, J. B.; Kästner, J. Adaptive Integration Grids in Instanton Theory Improve the Numerical Accuracy at Low Temperature. J. Chem. Phys. 2011, 134, 184107.
Figure 6. Arrhenius plot of the rate constant vs inverse temperature.
11.5 ms at a temperature of 50 K. The classical rate constant (which includes zero point vibrational energy, but no tunneling contributions) is many orders of magnitude smaller than the tunneling rate constant. The remaining slight temperature dependence of the tunneling rate constant shows that, while tunneling leads to synchronous concerted proton transfers, it is still aided by the finite thermal motion, a process known as temperature-assisted tunneling. The physical background is that in the reactant state, all molecules are further apart from each other than energetically favorable for the proton transfers. The contraction of the complex requires hardly any energy, but significant movements of atoms. This process is more efficient in classical movement than in tunneling of atoms down to very low temperature. The effect shows in long flat tails in the energy curves of the instantons at low temperature to the left in Figure 4 and in the still noticeable decrease of the tunneling rate at very low temperature.
■
CONCLUSION We have shown that quantum tunneling of atoms synchronizes the proton movements along Grotthuss chains. While the chosen example is a special case of a strongly exothermic reaction with an aluminum atom at the end of the transport chain, a similar mechanism may be relevant in liquid water. We have used instanton theory to calculate tunneling rate constants. Extended instanton solutions and corresponding tunneling rates were obtained in analytic cases and in a molecular example even above the crossover temperature Tc. While Tc is still a helpful concept to judge without expensive calculations if tunneling might play a role in a particular system, the present study shows that such analysis has to be interpreted with care.
■
REFERENCES
AUTHOR INFORMATION
Corresponding Author
*E-mail: [email protected]. Phone: +49 711 685 64473. Fax: +49 711 685 64442. Notes
The authors declare no competing financial interest.
■
ACKNOWLEDGMENTS J.R.F. and S.A.B. acknowledge the financial support of the University of Vigo (Contrato-Programa, 2012). S.A.B. acknowledges a F.P.U. grant from the Spanish Ministry of Education. J.K. would like to thank the German Research Foundation (DFG) for financial support of the project within the Cluster of 81
dx.doi.org/10.1021/jp411189m | J. Phys. Chem. A 2014, 118, 78−82
The Journal of Physical Chemistry A
Article
(22) de Grotthuss, C. J. T. Sur la décomposition de l’eau et des corps qu’elle tient en dissolution à l’aide de l’électricité galvanique. Ann. Chim. (Paris) 1806, 58, 54. (23) Hückel, E. Theory of the mobilities of the hydrogen and hydroxyl ions in aqueous solution. Z. Electrochem. 1928, 34, 546. (24) Bernal, J. D.; Fowler, R. H. A Theory of Water and Ionic Solution, with Particular Reference to Hydrogen and Hydroxyl Ions. J. Chem. Phys. 1933, 1, 515. (25) Agmon, N. The Grotthuss mechanism. Chem. Phys. Lett. 1995, 244, 456. (26) Cukierman, S. Et tu, Grotthuss! and other unfinished stories. Biochim. Biophys. Acta 2006, 1757, 876. (27) Knight, C.; Voth, G. A. The Curious Case of the Hydrated Proton. Acc. Chem. Res. 2012, 45, 101. (28) Kästner, J.; Carr, J. M.; Keal, T. W.; Thiel, W.; Wander, A.; Sherwood, P. DL-FIND: an Open-Source Geometry Optimizer for Atomistic Simulations. J. Phys. Chem. A 2009, 113, 11856. (29) Sherwood, P.; et al. QUASI: A general purpose implementation of the QM/MM approach and its application to problems in catalysis. J. Mol. Struct. (THEOCHEM) 2003, 632, 1. (30) Metz, S.; Kästner, J.; Sokol, A. A.; Keal, T. W.; Sherwood, P. ChemShell−a modular software package for QM/MM simulations. Comput. Mol. Sci. 2013, DOI: 10.1002/wcms.1163. (31) ChemShell, a Computational Chemistry Shell, see www. chemshell.org. (32) Zhao, Y.; Truhlar, D. G. The M06 suite of density functionals for main group thermochemistry, thermochemical kinetics, noncovalent interactions, excited states, and transition elements: two new functionals and systematic testing of four M06-class functionals and 12 other functionals. Theor. Chem. Acc. 2008, 120, 215. (33) Clark, T.; Chandrasekhar, J.; Spitznagel, G. W.; Schleyer, P. v. R. Efficient diffuse function-augmented basis sets for anion calculations. III. The 3-21+G basis set for first-row elements, Li-F. J. Comput. Chem. 1983, 4, 294. (34) Krishnan, R.; Binkley, J. S.; Seeger, R.; Pople, J. A. Selfconsistent molecular orbital methods. XX. A basis set for correlated wave functions. J. Chem. Phys. 1980, 72, 650. (35) Frisch, M. J. et al. Gaussian 09, Revision C.01. (36) Á lvarez-Barcia, S.; Flores, J. R. The interaction of Al atoms with water molecules: A theoretical study. J. Chem. Phys. 2009, 131, 174307.
82
dx.doi.org/10.1021/jp411189m | J. Phys. Chem. A 2014, 118, 78−82
