Srinivasan S. Iyengar Department of Chemistry and Department of Physics, Indiana University

1. Quantum wavepacket ab initio molecular dynamics:A computational approach for quantum dynamics in large systems Srinivasan S. Iyengar Department of Chemistry and Department of Physics, Indiana University Group members contributing to this work: Jacek Jakowski (post-doc), Isaiah Sumner (PhD student), Xiaohu Li (PhD student), Virginia Teige (BS, first year student) Funding:

2. Lipoxygenase: enzyme Ion (proton) channels Predictive computations: a few (grand) challenges • Bio enzyme: Lipoxygenase: Fatty acid oxidation • Rate determining step: hydrogen abstraction from fatty acid • KIE (kH/kD)=81 • Deuterium only twice as heavy as Hydrogen • generally expect kH/kD = 3-8 ! • weak Temp. dependence of rate • Nuclear quantum effects are critical • Conduction across molecular wires • Is the wire moving? • Reactive over multiple sites • Polarization due to electronic factor • Polymer-electrolyte fuel cells • Dynamics & temperature effects

3. Chemical Dynamics of electron-nuclear systems • Our efforts: approach for simultaneous dynamics of electrons and nuclei in large systems: • accurate quantum dynamical treatment of a few nuclei, • bulk of nuclei: treated classically to allow study of large (enzymes, for example) systems. • Electronic structure simultaneously described: evolves with nuclei • Spectroscopic study of small ionic clusters: including nuclear quantum effects • Proton tunneling in biological enzymes: ongoing effort

4. Hydrogen tunneling in Soybean Lipoxygenase 1: Introduce Quantum Wavepacket Ab Initio Molecular Dynamics Catalyzes oxidation of unsaturated fat Expt Observations • Rate determining step: hydrogen abstraction from fatty acid • Weak temperature dependence of k • kH/kD = 81 • Deuterium only twice as heavy as Hydrogen, • generally expect kH/kD = 3-8. • Remarkable deviation “Quantum” nuclei The electrons and the “other” classical nuclei

5. The “Quantum” nuclei The electrons and the “other” classical nuclei Quantum Wavepacket Ab Initio Molecular Dynamics [Distributed Approximating Functional (DAF) approximation to free propagator] Ab Initio Molecular Dynamics (AIMD) using: Atom-centered Density Matrix Propagation (ADMP) OR Born-Oppenheimer Molecular Dynamics (BOMD) References… S. S. Iyengar and J. Jakowski, J. Chem. Phys. 122 , 114105 (2005). Iyengar, TCA, In Press. J. Jakowski, I. Sumner and S. S. Iyengar, JCTC, In Press (Preprints: author’s website.)

6. 1. DAF quantum dynamical propagation Quantum Dynamics subsystem: Quantum Evolution: Linear combination of Hermite functions: The “Distributed Approximating Functional” is a banded, Toeplitz matrix Time-evolution: vibrationally non-adiabatic!! (Dynamics is not stuck to the ground vibrational state of the quantum particle.) Linear computational scaling with grid basis

7. Quantum dynamically averaged ab Initio Molecular Dynamics • Averaged BOMD: Kohn Sham DFT for electrons, classical nucl. Propagation • Approximate TISE for electrons • Computationally expensive. • Quantum averaged ADMP: • Classical dynamics of {RC, P}, through an adjustment of time-scales acceleration of density matrix,P “Fictitious” mass tensor ofP Force onP • V(RC,P,RQM;t) : the potential that quantum wavepacket experiences Ref.. Schlegel et al. JCP, 114, 9758 (2001). Iyengar, et. al. JCP, 115,10291 (2001).

8. Quantum Wavepacket Ab Initio Molecular Dynamics: The pieces of the puzzle The “Quantum” nuclei [Distributed Approximating Functional (DAF) approximation to free propagator] Simultaneous dynamics Ab Initio Molecular Dynamics (AIMD) using: Atom-centered Density Matrix Propagation (ADMP) OR Born-Oppenheimer Molecular Dynamics (BOMD) The electrons and the “other” classical nuclei S. S. Iyengar and J. Jakowski, J. Chem. Phys. 122 , 114105 (2005) J. Jakowski, I. Sumner, S. S. Iyengar, J. Chem. Theory and Comp. In Press

9. So, How does it all work? • A simple illustrative example: dynamics of ClHCl- • Chloride ions: AIMD • Shared proton: DAF wavepacket propagation • Electrons: B3LYP/6-311+G** As Cl- ions move, the potentialexperienced by the “quantum” proton changes dramatically. The proton wavepacket splits and simply goes crazy!

10. Spectroscopic Properties • The time-correlation function formalism plays a central role in non-equilibrium statistical mechanics. • When A and B are equivalent expressions, eq. (18) is an autocorrelation function. • The Fourier Transform of the velocity autocorrelation function represents the vibrational density of states.

11. Vibrational spectra including quantum dynamical effects • ClHCl- system: large quantum effects from the proton • Simple classical treatment of the proton: • Geometry optimization and frequency calculations: Large errors • Dimensionality of the proton is also important: • 1D, 2D and 3D treatment of the quntum proton provides different results. • McCoy, Gerber, Ratner, Kawaguchi, Neumark … • In our case: Use the wavepacket flux and classical nuclear velocities to obtain the vibrational spectra directly: • Includes quantum dynamical effects, temperature effects (through motion of classical nuclei) and electronic effects (DFT). In good agreement with Kawaguchi’s IR spectra References… J. Jakowski, I. Sumner and S. S. Iyengar, JCTC, In Press (Preprints: Iyengar Group website.)

12. The Main Bottleneck: The work around: Time-dependentDeterministic Sampling (TDDS) • Consider the phenol amine system Need the quantum mechanical Energy at all these grid points!! • However, some regions are more important than others? • Addressed through TDDS, “on-the-fly”

13. The Main Bottleneck: The quantum interaction potential The Interaction Potential: A major computational bottleneck The potential for wavepacket propagation is required at every grid point!! And the gradients are also required at these grid points!! Expensive from an electronic structure perspective • Quantum Dynamics subsystem: • AIMD subsystem (ADMP for example)

14. Time-dependent deterministic sampling 1)Importance of each grid point (RQM) based on: - large wavepacket density - r -potential is low- V • gradient of potential is high- V 2) So, the sampling function is: Ir, IV , IV’ --- adjust importance of each component

15. TDDS - Haar wavelet decomposition

16. Generalization to multidimensions - Haar wavelet decomposition

17. TDDS/Haar: How well does it work? The error, when the potential is evaluated only on a fraction of the points is really negligble!!! 1 mEh = 0.0006 kcal/mol = 2.7 * 10-5 eV Hence, PADDIS reproduces the energy: Computational gain three orders of magnitude!!

18. TDDS/Haar: Reproduces vibrational properties? The error in the vibrational spectrum: negligible These spectra include quantum dynamical effects of proton along with electronic effects!

19. Hydrogen tunneling in biological enzymes: The case for Soybean Lipoxygenase 1 Lipoxygenase: enzyme • Enzyme active site shown • Catalyzes the oxidation of unsaturated fat! • Rate determining step: hydrogen abstraction • Weak temperature dependence of k • Hydrogen to deuterium KIE is 81 • Deuterium is only twice as larger as Hydrogen, • Generally expect kH/kD = 3-8.

20. Lipoxygenase: enzyme Soybean Lipoxygenase 1: • A slow time-scale process for AIMD • Improved computational treatment through “forced” ADMP. • The idea is the donor atom is “pulled” slowly along the reaction coordinate • Bottomline: Donor acceptor distance is not constant during the hydrogen transfer process. • The donor-acceptor motion reduces barrier height

21. Soybean Lipoxygenase 1: Proton nuclear “orbitals”: Look for the “p” and “d” type functions!! s-type p-type d-type p-type These states are all within 10 kcal/mol Eigenstates obtained from Arnoldi iterative procedure

22. Reactant Transition State Eigenstates obtained using: • Instantaneous electronic structure (DFT: B3LYP) • finite difference approximation to the proton Hamiltonian. • Arnoldi iterative diagonalization of the resultant large (million by million) eigenvalue problem. For Deuterium, the excited proton state contributions are about 10% For hydrogen the excited state contribution is about 3% Significant in an Marcus type setting.

23. Transition state classical quantum H D

24. Conclusions and Outlook • Quantum Wavepacket ab initio molecular dynamics: Seems Robust and Powerful • Quantum dynamics: efficient with DAF • Vibrational non-adiabaticity for free • AIMD efficient through ADMP or BOMD • Potential is determined on-the-fly! • Importance sampling extends the power of the approach • In Progress: • QM/MM generalizations: Enzymes • generalizations to higher dimensions and more quantum particles: Condensed phase • Extended systems (Quantum Dynamical PBC): Fuel cells

25. Additional slides

26. Optimization of w‘(RQM) with respect to a,b,g a=1 (IY) b=3 (IYP) g=1 (IChi) RMS error of intrepolation during a dynamics within mikrohartrees

27. Computational advantages to DAF quantum propagation scheme • Free Propagator: is a banded, Toeplitz matrix: • Time-evolution: vibrationally non-adiabatic!! (Dynamics is not stuck to the ground vibrational state of the quantum particle.)

28. Quantum Wavepacket Ab Initio Molecular Dynamics: Working Equations Quantum Dynamics subsystem: Trotter Coordinate representation: • The action of the free propagator on a Gaussian: exactly known • Expand the wavepacket as a linear combination of Hermite Functions • Hermite Functions are derivatives of Gaussians • Therefore, the action of free propagator on the Hermite can be obtained in closed form: • Coordinate representation for the free propagator. Known as the Distributed Approximating Functional (DAF) [Hoffman and Kouri, c.a. 1992] • Wavepacket propagation on a grid

29. Spreading transformation We want to do potential evaluation for η fraction of grid • Density from ω(x) may be larger than current grid density- exceeding density is spread over low density grid area • for η 1 weightingω(x) should tend to 1 Grid point for potential evaluation are deteminned by integrating [N*w(x)] Interpolation of potential Version of cubic spline interpolation • based on on potentials and gradients • easy to generalize in multidimensions • general flexible form

30. Another example: Proton transfer in the phenol amine system Shared proton: DAF wavepacket propagation All other atoms: ADMP Electrons: B3LYP/6-31+G** C-C bond oscilates in phase with wavepacket Wavepacket amplitude near amine Scattering probability: References… S. S. Iyengar and J. Jakowski, J. Chem. Phys. 122 , 114105 (2005).

31. Potential Adapted Dynamically Driven Importance Sampling (PADDIS) : basic ideas The following regions of the potential energy surface are important: • Regions with lower values of potential • That’s probably where the WP likes to be • Regions with large gradients of potential • Tunneling may be important here • Regions with large wavepacket density Consequently, the PADDIS function is: The parameters provide flexibility