Chapter 7
Molecular Reaction Modeling from Ab-initio Molecular Dynamics 1
2
Peter E . Blöchl , Hans Martin Senn , and Antonio Togni
2
1
IBM Research Division, Zurich Research Laboratory, CH-8803 Rüschlikon, Switzerland Swiss Federal Institute of Technology ΕTH, Laboratory of Inorganic Chemistry, CH-8092 Zürich, Switzerland
Downloaded by NORTH CAROLINA STATE UNIV on August 7, 2012 | http://pubs.acs.org Publication Date: April 8, 1999 | doi: 10.1021/bk-1999-0721.ch007
2
A tutorial-like introduction to the methodology of transition state search within the framework of AIMD is given. We describe how to locate transition states, how to obtain dynamical reaction paths, and we discuss free energy integration. As an illustrating example we present results for the catalytic hydroamination of alkenes using {NiCl(PH ) } complexes.
+
3
2
Ab-initio molecular dynamics (AIMD), invented 1985 by Car and Parrinello (/), has established itself as a useful tool for the study of catalytic processes. The methodol ogy used in connection with AIMD takes advantage of a number of concepts which originated in solid state physics and adapts others from chemistry. While being based on state-of-the-art density functional methodology (2, 3), AIMD opens a wide range of options that have been accessible in the past only in combination with simpler semi-empirical or molecular mechanics methods. It thus expands the spectrum of theoretical chemistry into new directions. The calculations presented here are based on the projector-augmented wave (PAW) method (4), which has been developed to make the full wave functions avail able within the AIMD approach, avoiding the usual pseudopotential approximations (5, 6). Quantities depending on the charge density near the nucleus, such as electric field gradients, are therefore directly accessible with high accuracy (Petrilli, H. M.; Blochl, P. E.; Blaha, P.; Schwarz, K. Phys. Rev. Β 1998, in press). We describe in a tutorial-like way aspects of the methodology related to transi tion states as used within the framework of AIMD. We discuss how to locate transi tion states, how to analyze the dynamics of the reaction event, and how to perform free energy integrations. These methods have been used earlier in connection with molecular mechanics, semiempirical, and AIMD simulations, and they belong now, in variants, to the toolset of state-of-the-art AIMD simulations performed by ourselves and others.
88
© 1999 American Chemical Society
In Transition State Modeling for Catalysis; Truhlar, D., et al.; ACS Symposium Series; American Chemical Society: Washington, DC, 1999.
89
Downloaded by NORTH CAROLINA STATE UNIV on August 7, 2012 | http://pubs.acs.org Publication Date: April 8, 1999 | doi: 10.1021/bk-1999-0721.ch007
Methodology A l l our calculations are based on density functional theory (2, 3). In particular, the results presented here are obtained using the gradient corrections due to Becke (7) and Perdew (8). We use the general framework of A I M D (7), termed the "Car-Parrinello method" after its inventors. The approach used here is based on the fictitious Lagrangian methodology (7, 9), in which the electron wave functions are propagated in time via a differential equation. If properly done, the wave functions follow the nuclei adiabatically on, or in practice very close to, the Born-Oppenheimer surface, that is with the electron wave functions in the instantaneous ground state. Both atoms and wave functions follow dynamical equations of motion: (1)
Here, E is the Kohn-Sham total energy functional, which depends on the atomic po sitions R; and on the one-particle electron wave functions | ^ ) ; M , are the nuclear masses, and is the fictitious mass of the wave functions. H is the one-particle Hamiltonian acting on the wave functions, and A are the Lagrange parameters for the constraint that all wave functions are orthonormal. The analogy between the two differential equations becomes evident, when we express the Hamiltonian via a partial derivative H\%) = 3 £ / 3 ( ^ J of the Kohn-Sham total energy functional. The advantage of this approach lies in the absence of self-consistency cycles as the atoms move. Other variants of AIMD (70) use self-consistency cycles, while still exploiting the information of the previous steps. Those approaches allow larger time steps for the atom dynamics, but at the cost of self-consistency cycles. The A I M D approach can be implemented in various electronic structure meth ods. Most commonly, pseudopotentials are used (5, 6). The electronic structure prob lem in our work, on the other hand, has been solved using the PAW method (4). It is an all-electron method in the sense that the full, valence and core, wave functions and densities are constructed to evaluate the total energy and forces. However, we use the frozen core approximation (77), that is, we import the core densities of an isolated atom into the molecular calculation. This turns out to be a mi nor approximation, if the exact core densities are used, because then the errors in the total energy are only of second order in the deviation of the true and the atomic core density. Note, however, that the frozen core states cannot be identified on a one-toone basis with one-particle states of the molecule. If one-particle states are needed, the consistent approach is to diagonalize the Hamiltonian, whose matrix elements are obtained from core and valence states resulting from the frozen core calculation. Such a transformation mixes individual states, but leaves the total energy, charge density, and forces untouched. The main characteristics of the PAW method are the choice of the basis functions and the procedure to calculate total energy and forces. We use basis functions com posed of plane waves, which are augmented by atomic "orbitals" near the nucleus in order to incorporate the proper nodal structure of the wave functions. A wave func tion takes the form n
mn
(2) RJ,m,q
In Transition State Modeling for Catalysis; Truhlar, D., et al.; ACS Symposium Series; American Chemical Society: Washington, DC, 1999.
90 where
are the so-called pseudo wave functions, which are the variational pa
rameters, that is they contain the expansion coefficients c
Gn
as ^„(r) = ^ e G
l G r
CG,„.
The true wave function is obtained by adding and subtracting one-center expansions in terms of partial waves \R i, , ) and |^,/, ,^), which are atom-centered radial %
m q
m
Downloaded by NORTH CAROLINA STATE UNIV on August 7, 2012 | http://pubs.acs.org Publication Date: April 8, 1999 | doi: 10.1021/bk-1999-0721.ch007
functions multiplied by spherical harmonics. The index R denotes a nuclear site, /, m are angular momentum indices, and q is an additional counter. In the following, we will collect the indices /, m, q into one integer counter i or j. The all-electron partial waves |