Approximate Quantum Dynamics using Ab Initio ... - ACS Publications

Subscriber access provided by Fudan University

Article

Approximate Quantum Dynamics using Ab Initio Classical Separable Potentials: Spectroscopic Applications Barak Hirshberg, Lior Sagiv, and Robert Benny Gerber J. Chem. Theory Comput., Just Accepted Manuscript • Publication Date (Web): 08 Feb 2017 Downloaded from http://pubs.acs.org on February 8, 2017

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a free service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are accessible to all readers and citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

Journal of Chemical Theory and Computation is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 33

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

Approximate Quantum Dynamics using Ab Initio Classical Separable Potentials: Spectroscopic Applications Barak Hirshberg1, Lior Sagiv1 and R. Benny Gerber1,2* 1) Institute of Chemistry and the Fritz Haber Center for Molecular Dynamics, the Hebrew University, Jerusalem 9190401, Israel. 2) Department of Chemistry, University of California, Irvine, California 92697, USA KEYWORDS quantum dynamics, classical separable potentials, AICSP, vibrational spectroscopy, ab initio potential energy surfaces, amino acids, guanine-cytosine base pair

ABSTRACT Algorithms for quantum molecular dynamics simulations which use directly ab initio methods have many potential applications. In this paper, the Ab Initio Classical Separable Potentials method (AICSP) is proposed as the basis for approximate algorithms of this type. The AICSP method assumes separability of the total time-dependent wavefunction of the nuclei, and employs mean-field potentials that govern the dynamics of each degree of freedom. In the proposed approach, the mean-field potentials are determined by classical ab initio molecular dynamics (AIMD) simulations. The nuclear wavefunction can thus be propagated in time using the effective potentials generated "on the fly". As a test of the method for realistic systems,

ACS Paragon Plus Environment

1


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 2 of 33

calculations of the stationary anharmonic frequencies of hydrogen stretching modes were carried out for several polyatomic systems, including three amino acids and the guanine-cytosine (G-C) pair of nucleobases. Good agreement with experiments was found. The method scales very favorably with the number of vibrational modes, and should be applicable for very large molecules, e.g. peptides. The method should also be applicable for properties such as vibrational linewidths and lineshapes. Work in these directions is underway.

A. Introduction Ab initio molecular dynamics (AIMD) simulations are a powerful tool for exploring properties of molecular systems, e.g., structure and bonding of new materials1, chemical reactions in gas and in condensed phases2, studying biomolecular processes3-5 etc. These methods include the direct dynamics or Born-Oppenheimer molecular dynamics (BOMD) approach6, in which the electronic time independent Schrodinger equation (TISE) is solved to obtain the forces acting on the nuclei at every time step, as well as the celebrated Car-Parrinello molecular dynamics (CPMD) method7. In these simulations, the motion of the nuclei is described classically by solving Newton's equations of motion using appropriately sampled initial conditions. The classical treatment is an excellent approximation for heavy atoms and relatively high temperatures. However, at low temperatures, or for reactions which involve the motion of light atoms, nuclear quantum effects (zero-point energy, tunneling, interference etc.) may become significant. Ideally, one would like to solve the full time-dependent Schrodinger equation (TDSE) for the nuclei, using a suitable Born-Oppenheimer (BO) potential energy surface (PES). However, this task is extremely difficult for systems having more than a few atoms due to the exponential scaling of the computational complexity with the number of degrees of freedom


2

Page 3 of 33

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


involved8. Several important methods were developed in order to reduce the computational cost of numerical solutions to the TDSE8-13. However, despite these important contributions, state-ofthe-art numerical solution of the TDSE in full dimensionality is still limited to systems of several atoms only. As a result, there are many potential applications for a quantum mechanical method which: 1) Uses ab initio PES, as is often essential for accuracy. Direct use of PES, if possible is preferable. 2) Is a time-dependent method, which allows the treatment of non-equilibrium processes. 3) Has the capacity to be applied to systems having a large number of degrees of freedom. This paper will present such a method, the ab initio Classical Separable Potentials method (denoted AICSP). This is an approximate but accurate time-dependent method, which is applicable to large systems and uses first-principles potentials directly. Several excellent methods are in use today which possess some, but not all, of the properties described above. An elegant review on time-dependent quantum dynamics methods is by Makri14. Here, we discuss briefly some of the methods, and highlight their strengths and limitations. Time-Dependent Self-Consistent Field (TDSCF) and CSP. Possible alternatives to the exact numerical solution of the nuclear TDSE include mean-field theories14-15. Perhaps the first meanfield theory was the time-dependent self-consistent field method (TDSCF), often referred to as the time-dependent Hartree (TDH) method16. In this method, the total nuclear wavefunction is assumed to be completely separable and a set of coupled equations is obtained in which each mode evolves under an average one-dimensional potential. This potential is the quantummechanical average of the full PES over all other nuclear degrees of freedom, and thus, a selfconsistent solution is needed (see Section B for details). The origins of the approach lie in the very early days of quantum mechanics, with first applications to chemical dynamics dating from 198216 for the dissociation of van der Waals clusters. Semiclassical and classical variants of


3


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 4 of 33

TDSCF were also introduced16-17. Since its development, the method has been applied, in good agreement with numerically exact solutions, to several problems such as photodissociation processes18-21, unimolecular decomposition22, vibrational predissociation23 and others24-26. Despite the significant computational advantage over numerical solution of the full TDSE, TDSCF is also limited to rather small systems. In addition, this method was not applied using ab initio potentials. A closely related approach, which proved very useful in describing larger systems is the Classical Separable Potentials (CSP) method27. Using this method, the singlemode potentials are not obtained by quantum-mechanically averaging the full PES over all other modes of vibration but rather by an average obtained from classical trajectories simulations (see Section B for details). This method was applied, with encouraging results, to ultrafast dynamics in large noble-gas clusters28-30, electron photodetachment from C60- 31 and more32-33. Despite the additional approximation involved in the CSP method, it was shown that for short timescales quantitative agreement between CSP and TDSCF is obtained33. Almost all of CSP applications were done using empirical force fields. One exception is the application by Knospe and Jungwirth31 in which the electron photodetachment from C60- was studied using a PES calculated "on the fly" using a DFT method. As mentioned previously, it is the goal of this paper to extend the CSP method to use ab initio PES. Multi-Configuration Time-Dependent Hartree (MCTDH). Although TDSCF and CSP have proved to be accurate for short time dynamics of a variety of processes, in some cases the single configuration ansatz used in the TDSCF and CSP methods is just not sufficient34. The inclusion of multiple configurations in the wavefunction can be typically done by the multi-configurational time-dependent Hartree (MCTDH) method35-38. Although this method is more efficient than exact numerical solution of the TDSE, it is still computationally much more demanding than


4

Page 5 of 33

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


TDSCF. For small systems of several atoms, MCTDH can effectively provide rigorous quantum solutions while for large systems, it is a high-level approximation. Very efficient algorithms in recent years significantly reduce the computational cost39-41 and allow applications of this method to increasingly large systems. However, most applications are done using empirical potentials or model Hamiltonians. Few applications were also done with analytic potentials fitted to ab initio calculations for small molecules42-43. Recently, the direct dynamics variational multiconfigurational Gaussian technique (DD-vMCG) was suggested44 which is based on the MCTDH method for quantum nuclear dynamics and uses ab initio potentials. So far, the method has been applied to systems with 10-20 atoms45-46. Path Integral methods. Other methods which are very useful in incorporating quantum mechanical effects for the nuclei are centroid molecular dynamics (CMD) and ring-polymer molecular dynamics (RPMD)14,

47-49

, which are based on the path-integral formulation of

quantum mechanics. An excellent review was published recently by Habershon et al.47 which discusses the theory behind these methods, their advantages and limitations. Already several years ago, these methods were used directly with ab initio PES, but for relatively small systems

such as hydrazine, and 50-53. Very recently, efficient schemes for ab initio path-integral

molecular dynamics were suggested54-56. These even allow the simulation of liquid water using

RPMD with ab initio potentials at significantly reduced computational cost54. While these algorithms seem very promising, most applications of ab initio RPMD and CMD methods are to calculate equilibrium properties, such as structure and energetics. We note that these methods were also extended to provide approximate time-correlation functions47, 49, 57 which can be used to calculate IR absorption spectra, diffusion coefficients, reaction rates and more. However, it is not straight-forward, for example, to apply these methods for dynamically simulating chemical


5


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 6 of 33

reactions. In addition, while the approximate time-correlation functions obtained from CMD or RPMD yield good agreement for diffusion coefficients and reaction rates, the methods are less successful in calculating absorption spectra of high-frequency vibrational modes47,

58-59

.

Alternative methods which are able to both simulate chemical reactions dynamically in time as well as calculate the vibrational absorption spectra in the high-frequency range may thus be very useful. As a result, this paper will focus on high-frequency hydrogen stretching vibrational modes. We limit the scope of this paper to quantum dynamical methods. Semiclassical methods such as the Semiclassical Initial Value Representation method (SC-IVR)60-61 or the Gaussian wavepacket method17, 62-63 can in principle treat systems with a large number of degrees of freedom. Several recent developments applied these methods using ab initio potentials64-66. While most applications were done for small systems such as CO2 and formaldehyde, Vanicek et al. calculated66 the emission spectra of pentathiophene with 105 vibrational degrees of freedom. To summarize, there is currently no method for quantum dynamical calculations which can be routinely applied to large systems and use directly ab initio PES. The main goal of this paper is to present a new method which combines the CSP method with ab initio PES obtained from quantum chemistry. We will present results using both Möller-Plesset second order perturbation theory (MP2)67 potentials and a long-range corrected hybrid density functional with an empirical dispersion correction, B97X-D68.

In applying the new method, testing is of the essence. Stationary vibrational spectroscopy presents a very useful way of doing so. While a time-dependent method can treat many nonstationary processes, few of these make a strict comparison as stationary vibrational


6

Page 7 of 33

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


spectroscopy, since comparison with extremely accurate experimental results and high-level computational methods is possible. The point of view taken here is that if the method can compete with computational methods for vibrational spectroscopy in accuracy and systems size when comparing to experimental results, it strengthens the assumption that it will also be useful for future, time-dependent applications. We will test AICSP for anharmonic vibrational spectroscopy of amino acids (glycine, alanine and proline), the guanine-cytosine (G-C) base pair and HNO3. These systems are on their own not trivial when ab initio potentials are used (for example, proline has 45 vibrational modes) and can also serve as a benchmark before an attempt to tackle very large systems such as long peptides. So far, neither CSP nor TDSCF were applied to stationary vibrational spectroscopy, which is a quantitatively demanding test. The paper is organized as follows: In Section B, the theory of TDSCF and CSP methods is outlined. Section C reports the details of the computational implementation of the theory and Section D reports the results for amino acids, G-C base pair and HNO3. The last Section summarizes and presents the conclusions of the work. B. Theory: The TDSCF method and its relation to the CSP method Within the TDSCF method, the total nuclear wavefunction is assumed to be separable, typically in normal coordinates (1)

Χ , … , , = ,

Where Χ is the total nuclear wavefunction, are the single-mode wavefunctions, is the number of vibrational modes of the system and are the mass-weighted normal modes of ACS Paragon Plus Environment

7


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 8 of 33

vibration. Using this ansatz, and applying the Dirac-Frenkel variational principle69-70, the TDSCF equations are obtained: (2)

−

ℏ , + 2

!"

, , = #ℏ

, ; & = 1, … ,

Where , are the same as the single mode functions , up to a purely timedependent phase factor16.

!"

, is the average potential for mode , which is given by

the quantum mechanical average of the full PES, , … , , over all other modes of vibration. In this paper, , … , is the ab initio PES used.

(3)

!"

, = ( ) * ) , +- , … , - ) * ) , +. ),

),

As a result, a self-consistent solution to the TDSCF equations is needed. The TDSF equations strictly conserve the total energy and the norm of the total nuclear wavefunction16. It was found previously, in test cases, that the TDSCF approximations give good agreement with numerically exact quantum solutions16, 22-24, especially for short timescales. Calculating the average potential requires multidimensional integration and is the most significant computational bottleneck of the TDSCF method. This becomes very difficult to achieve, even for systems of moderate size ( ≈ 10, unless pairwise coupling between the modes is assumed. The CSP method can be thought of as introducing a further approximation to TDSCF, in which the multidimensional

integrals are approximated by an average value obtained from classical trajectories simulations. The classical averaging must be performed over a set of trajectories which describe the classical analogue for the quantum mechanical calculation. For example, if one is interested in the energy


8

Page 9 of 33

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


transfer from one vibrational mode to another, the initial conditions should be sampled

accordingly. More specifically, using 123) classical trajectories, given by 4 5 65 ,…, 789: for every normal mode , the average single-mode potentials in the CSP method are obtained by

(4)

(5)

5 , = 5 , … , , … , 5 ! ; , =

1

123)

789:

< 5 ,

5

where is any point on the & normal coordinate. Once the average single-mode potentials are

obtained, the time-dependent uncoupled equations can be solved for each mode separately according to the following equations (6)

ℏ , , ! ; , , = #ℏ ; & = 1, … , − + 2

We note that while the average single-mode potentials are now obtained in an entirely classical manner, all modes of vibration are completely coupled to one another in determining the effective potential. This removes the necessity of solving the equations self-consistently, and significantly reduces the computational cost. In addition, due to the time-dependence of the potential, energy can flow between different modes. While the non-self-consistency may in principle affect the accuracy, it was shown in test calculations that the CSP method gives results with very good agreement with TDSCF, for systems which exhibit moderate quantum effects and for short timescales33. In addition, the approximation that the single-mode effective potentials can be evaluated using classical trajectories was examined previously and the potentials obtained in this way were in very agreement with effective single-mode potentials obtained within the


9


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 10 of 33

TDSCF framework33. Unlike TDSCF, the CSP equations do not strictly conserve the total energy. However, it was shown previously that for short timescales the change in total energy is small and does not seem to introduce significant errors33. In this paper we test AICSP by calculation of stationary vibrational frequencies. Linewidths or lineshapes will be addressed in the future. The key issue for calculating the stationary spectra from the CSP approach is to determine the appropriate effective potentials. Our aim is to describe the system for zero temperature, where quantum mechanically the system has zero point energy for the vibrational ground state. Choosing zero vibrational energy in describing the classical analogue is inadequate since it takes all modes, except the one being considered for frequency calculations, as being frozen at equilibrium positions. All the effects of dynamical couplings between different vibrations will then be missing. Introducing an energy to represent the ZPE effect in the effective potential calculations seems a better approach. We chose to give equal energy to each mode, in the spirit of the classical equipartition behavior. Further, we choose the harmonic quantum ZPE value of the lowest frequency mode in this role. This is to avoid large energy flow effects between different modes, that do not take place physically at T=0 K. Using a much larger value for the classical analogue of ZPE may also lead to the unphysical behavior of conformational changes at zero temperature. After using this procedure for initial conditions sampling, which approximately describe a stationary state of the system, a reasonable approach to obtain vibrational eigenstates is to average the time-dependent single-mode

potentials ! ; , over some time period = to obtain a time-independent single mode

potential


10

Page 11 of 33

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


(7)

?! ; = V

1 B !; @ , A = C

Then, it is trivial to numerically solve the one dimensional TISE for every mode and obtain anharmonic wavefunctions and frequencies. (8)

−

ℏ D ?! ; D = E D +V 2

In the future, we intend to extend the method and obtain vibrational lineshapes for similar systems by solving the time-dependent single-mode equations. However, stationary vibrational spectroscopy is a useful test for the purpose of evaluating the suitability of key components of the AICSP method, such as using effective potentials obtained from classical trajectories. We stress here that both TDSCF and CSP are based on the assumption that the total nuclear wavefunction can be approximated by a separable ansatz. In this implementation of the AICSP method, normal coordinates are used. As a result, the method is currently not adequate for extremely anharmonic systems, for which separability in normal modes is not a good approximation, or for cases when the dynamics in time are not limited to the vicinity of a single minimum on the PES. It is possible that these limitations may be removed by a better choice of coordinates, such as internal or local coordinates, which were successfully used in Vibrational Self-Consistent Field71 (VSCF) calculations72-73. C. Computational Details


11


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 12 of 33

Prior to the CSP simulations a geometry optimization is performed, followed by a calculation to determine the harmonic vibrational frequencies and normal modes to verify that the structure is a true minimum on the PES with no imaginary frequencies. In addition, the harmonic frequencies and normal modes are used to sample the initial conditions for the classical MD simulations, as will be explained shortly. All geometry optimizations, single point energy calculations and BOMD simulations were performed using Q-Chem 4.374 with either MP267 level of theory or B97X-D68 density functional theory, both using the aug-cc-pVDZ basis set75.

The implementation of the AICSP method is comprised of several different steps: 1. Sampling of initial conditions and performing classical MD simulations Initial values for each mode were randomly selected and the momentum for each normal mode was determined according to the single mode initial energy. As explained in section B, the initial mode energy was taken to be the ZPE of the lowest frequency mode of the system. For

simplicity, the initial values of were selected from the probability distribution of the classical

harmonic oscillator. Unphysical values of were also rejected based on the classical turning

points of the harmonic oscillator. The choice of the harmonic approximation at this point is not unique or essential. In any case, the effective potentials are calculated from the actual trajectories that include anharmonic effect in full. After sampling 123) initial conditions, classical BOMD simulations are performed using Q-

Chem 4.374 employing the Fock matrix extrapolation technique76. A time step of 0.48 fs was

used in the simulations, in order to accurately describe the light atoms motion. Each trajectory ran for a total time of 1.2 ps. We limit ourselves in this paper to the stretching frequencies of NH and O-H bonds, for which the typical vibrational period is approximately 10 fs. As a result,


12

Page 13 of 33

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


this total simulation time allows averaging over about 120 periods of the stretching mode. We stress that this method is not limited to high frequency stretching modes but a longer averaging time might be required for softer vibrational modes. In addition, for slower processes which involve soft vibrational modes, a larger number of classical trajectories may be required in order to improve the quality of the effective potentials. 2. Building the average single-mode potentials from the classical trajectories After 123) classical trajectories are obtained, we transform from Cartesian to normal

coordinates at each time step to obtain 4 5 65 ,…, 789: for every normal mode . Next, these values are used to build the single-mode average potential ! ; , on a spatial and

temporal grid. The spatial grid for was constructed from F2GH equally spaced points in the

range [−3K , 3K ], where K = MN is the standard deviation of the probability density for the ℏ

O

harmonic ground state of mode . The temporal grid was composed of every P1QR time steps of the classical trajectories.

In order to build the potential ! ; , for every mode according to Equations ( 4 )-( 5 ),

the following calculations were performed: For a given trajectory S, and every temporal grid

point, a set of F2GH ab initio single-point energy calculations was performed for every grid value of . In each energy calculation, all modes ), are fixed to their values from the classical simulations at the specific time step, for the specific trajectory, )5 . Together, these points comprise 5 , . This procedure is repeated for all trajectories and ! ; , is obtained as

the average of 5 , over all trajectories using equal weights. For the results presented in

this paper we used 123) = 10, F2GH = 17 and P1QR = 100. The numerical error due to ACS Paragon Plus Environment

13


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 14 of 33

changes in the values of these parameters is found to be smaller than ± 15 WXY . We also note

that in this approach, where the classical trajectories are performed as a preliminary step to building the potential, every single point energy calculation is entirely independent of the others. As a result, the potential building step can be trivially parallelized over a large number of CPUs. This make the possibility of using the CSP method for calculating O-H and N-H stretching frequencies in large systems more likely. However, the size of the systems treated will be limited by the ability to perform classical AIMD simulations of the system. 3. Averaging of the single-mode potentials over a time period Z and solution of TISE

Once ! ; , is obtained, it is averaged over a time period =, as shown in Equation ( 7 ).

Then, the one dimensional TISE given by Equation ( 8 ) is solved on a spatial grid using a second

order central finite difference approach for the kinetic energy. The eigenvalues and eigenstates

give the anharmonic frequencies and wavefunctions for each mode. = is determined so that the frequencies obtained are least sensitive to the averaging period. We found that for C-H, N-H and

O-H stretching modes in amino acids = ≈ 0.6 − 1 ]^ gives an averaging error of approximately

± 15 WXY .

D. Results and Analysis We employ the AICSP method to study the vibrational spectroscopy of amino acids and HNO3. Amino acids were chosen because they constitute the basic building blocks of proteins, and can serve as a reliable benchmark for the anharmonic method before treating much larger systems, such as large peptides or proteins. Vibrational spectroscopy of proteins is a powerful tool which can give information regarding the protonation states of side chains, strength of hydrogen bonds


14

Page 15 of 33

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


and more77-78. Extensive theoretical and computational results are available for these amino acids, which renders these systems a suitable framework for testing AICSP. HNO3 was chosen as a simple model system where accurate data exists for overtone vibrational spectroscopy. The different structures studied in this paper are presented in Figure 1. Each amino acid has several conformers. However, not all of them necessarily contribute to the vibrational spectrum under the experimental conditions. Since the experiments for the isolated amino acids are done in cold rare-gas matrices, where the temperature is low, only several low lying conformers are considered. Specifically, we present results for conformers IA and IIA of alanine and proline which differ mainly by the fact that, unlike the IA conformers, the IIA conformers have an internal hydrogen bond. Results are also presented for conformer I of glycine and for HNO3.


15


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 16 of 33

Figure 1. List of molecules considered in this study: HNO3, conformers IA and IIA of alanine (AlaIA and AlaIIA, respectively) and proline (ProIA and ProIIA, respectively) and conformer I of glycine (GlyI). Blue, red, black and white sphere represent nitrogen, oxygen, carbon and hydrogen atoms respectively. Proline. The gas-phase vibrational spectroscopy of proline has been studied previously experimentally79-80 and theoretically81. Mainly, Stepanian et al. have measured the IR spectra of proline at 14K in a solid argon matrix79-80 and identified two conformers in the matrix, ProIA and ProIIA. They have assigned the absorption peaks corresponding to each conformers based on frequencies calculated using the harmonic approximation at the B3LYP/aug-cc-pVDZ level of theory. As a result, we consider the same two conformers of proline: ProIA and ProIIA. Both conformers have one O-H and one N-H stretching modes, but they differ in the existence of an internal hydrogen bond between the hydrogen of the carboxylic acid and the nitrogen of the amine for the ProIIA conformer. As a result, the O-H stretching mode of the ProIIA is significantly red-shifted (by 534 cm-1) with respect to the ProIA O-H stretching frequency. In order to be able to determine different conformers of amino acids, for example, the computational method should reproduce this significant redshift. Figure 2 shows the results obtained for both proline conformers using different potentials and different anharmonic treatments. We compare the AICSP method to the harmonic approximation, as well as to the VSCF method71. Specifically, results from Brauer et al.81 using the VSCF procedure with hybrid PM3/MP2 potentials are used. In this procedure, the harmonic normal coordinates obtained using PM3 are scaled by the ratio between the MP2 and PM3 harmonic frequencies before being used in the calculation of VSCF anharmonic frequencies. This was shown to give improved agreement with experimental results. In addition, results for the AICSP method using both MP2 and


16

Page 17 of 33

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


B97X-D levels of theory are shown. It is clear that the agreement with experimental

frequencies for AICSP is far better than the harmonic approximation, which is not surprising, since it is an anharmonic method. Using MP2 potentials, the agreement for the N-H stretching modes of both conformers is very good with deviations that are smaller than 1%. This is also true for the frequencies calculated using B97X-D. For the O-H stretch of the ProIA conformer

excellent agreement is also seen for MP2 PES (deviation of 4 cm-1 or 0.1%). For the internally

hydrogen bonded O-H stretch, the deviation of MP2 from the experimental frequency is admittedly larger (approximately 4.1%). This may be due to the fact that hydrogen bonded structures are more sensitive to the interactions between modes. The CSP mode displacements possibly underestimate these. However, AICSP correctly results in a very significant red-shift which lowers the frequency for the ProIIA O-H stretch by 405 cm-1 with respect to the O-H stretch of ProIA. We note that MP2 shows somewhat better agreement than B97X-D for the O-

H stretches. Specifically, for the hydrogen bonded O-H mode, B97X-D shows a deviation of

6.3% from the experimental value. Generally, we see that the results of AICSP and VSCF are of similar quality. This is reasonable since VSCF is a time-independent version of the TDSCF method, upon which CSP is based.


17


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 18 of 33

Figure 2. Correlation plots comparing calculated and experimental frequencies (in cm-1) for O-H and N-H stretching modes for proline conformers. The black line represents perfect correlation between theory and experiment. Blue triangles represent the harmonic approximation using MP2 PES. Green circles represent VSCF results using hybrid PM3/MP2 potentials with the TZP basis set81. Black and red squares show results from AICSP theory using B97X-D and MP2

potentials, respectively.

Alanine and Glycine. Alanine and glycine were also studied experimentally by Stepanian et al.82-83 The spectra for glycine and alanine were collected at 13K and 14K, respectively, in an argon matrix. Similarly to proline, only two conformers of alanine were identified in the experiments using harmonic frequencies calculated at the B3LYP/aug-cc-pVDZ level of theory.


18

Page 19 of 33

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


As a result, two conformers of alanine were again considered: AlaIA and AlaIIA which also differ by the existence of an internal hydrogen bond for the AlaIIA conformer, similar to the one in ProIIA. For the AlaIIA isomer, the redshift due to the internal hydrogen bond is 367 cm-1 which is large, but somewhat smaller than for ProIIA. In the case of glycine three conformers were identified in the experiments; however, the major contribution was due to GlyI, which is the lowest energy conformer. In addition, conformer GlyIII was found to convert to GlyI upon matrix annealing. For glycine, only the GlyI conformer is considered, which does not have an internal hydrogen bond. Alanine and glycine were also studied theoretically using the VSCF approach with hybrid potentials, similarly to proline81. Figure 3 shows the results obtained for alanine and glycine using AICSP, in comparison to the harmonic approximation, as well as VSCF. For alanine, two potentials were again tested using the MP2 and B97X-D levels of

theory. For the MP2 potentials, excellent agreement is seen for the free O-H stretch of GlyI and AlaIA (deviations smaller than 6 cm-1 or 0.2%). For the hydrogen bonded O-H stretching mode of AlaIIA better agreement is seen, using MP2 potentials, than for ProIIA (deviation is 74 cm-1

or 2.3%) but still not as good as for the free O-H modes. However, also in this case a significant

red shift of 289 cm-1 is reproduced using MP2 potentials. B97X-D again shows somewhat

inferior agreement with experimental frequencies for the O-H vibrations. This is emphasized for

the O-H stretch of the AlaIIA where B97X-D deviates by more than 200 cm-1 from the

measured value and shows a much smaller redshift of 130 cm-1. Lastly, also for the case of alanine and glycine, VSCF and AICSP give results with similar deviations from experiments.


19


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 20 of 33

Figure 3. Correlation plots comparing calculated and experimental frequencies (in cm-1) for O-H stretching modes for alanine and glycine conformers. The black line represents perfect correlation between theory and experiment. Blue triangles represent the harmonic approximation using MP2 PES. Green circles represent VSCF results using hybrid PM3/MP2 potentials with the TZP basis set81. Black and red squares show results from AICSP theory using B97X-D and

MP2 potentials, respectively.

HNO3. The gas-phase vibrational spectroscopy of HNO3 was studied previously both experimentally and theoretically84-88. Here, we use it as a benchmark to test the accuracy of AICSP for overtones, as well as fundamentals. Figure 4 shows the comparison of frequencies calculated using the AICSP with PES obtained using MP2 level of theory, as well as B97X-D. ACS Paragon Plus Environment

20

Page 21 of 33

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


For HNO3, AICSP is compared to the VSCF-PT2 method89, which incorporates correlation corrections to VSCF using perturbation theory. Specifically, we compare to results obtained using PES at the MP2 level of theory88. In this case, excellent agreement for the fundamental

transition is seen for B97X-D, (deviation of 15 cm-1) while a deviation of almost 100 cm-1 is

seen for MP2. The situation is reversed for the first overtone where MP2 shows a very small

deviation of less than 1%. More interestingly, while a significant absolute error (approximately 440 cm-1 using MP2) is seen for the second overtone, it is still small on relative terms (4.3%). The agreement for VSCF-PT2 is excellent for the fundamental as well as both overtones. However, we note that using VSCF-PT2 to calculate overtone frequencies requires the same computational time as calculating a fundamental while using the method presented in this paper, the overtone frequency is trivially obtained from the solution to the single-mode TISE. In order to obtain reliable overtones using the AICSP method it is essential to make sure that the number of grid points is sufficient for describing the excited state wavefunction.


21


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 22 of 33

Figure 4. Correlation plots comparing calculated and experimental frequencies (in cm-1) for the

fundamental (0 → 1) and first (0 → 2) and second (0 → 3) overtones of the O-H stretch mode of

HNO3. The black line represents perfect correlation between theory and experiment. Blue triangles represent the harmonic approximation using MP2 PES. Green circles represent VSCFPT2 results using MP2 potentials with the TZP basis set88. Black and red squares show results from AICSP theory using B97X-D and MP2 potentials, respectively.

Proof of concept for large systems. In order to demostrate that the method is indeed useful for large systems, we calculated the stationary vibrational frequencies for C-H, N-H and O-H stretching modes of the guanine-cytosine (G-C) base pair which has 81 vibrational degrees of freedom. Specifically, the calculation is done for the isomer which is detected in the gas phase90


22

Page 23 of 33

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


in which the cytosine is in the enol form and there is a hydrogen in position 7 on guanine (see Figure 5). The results obtained using B97X-D potentials with a smaller basis set (cc-pVDZ75)

are presented in Figure 5.

Figure 5. Correlation plots comparing calculated and experimental frequencies (in cm-1) for the G-C base pair. The black line represents perfect correlation between theory and experiment. Blue triangles represent the harmonic approximation using B97X-D PES. Green circles represent ACS Paragon Plus Environment

23


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 24 of 33

VSCF-PT2 results using hybrid PM3/RI-MP2 potentials with the cc-pVDZ basis set91. Red

squares show results from AICSP theory using B97X-D potentials with the cc-pVDZ basis set. Panel (a) shows results for C-H and O-H stretching modes. Panel (b) shows results for N-H

stretching modes of the cytosine moiety. Panel (c) shows results for N-H stretching modes of the guanine moiety. Panel (d) shows the optimized structure for the enol form of the G-C complex. The results are compared to results obtained using the VSCF-PT2 method89 with hybrid PM3/RIMP2 potentials91, to results obtained for the harmonic approximation, and to experimental measurements of Abu-Riziq et al.90. Overall good agreement is seen for all modes with a mean absolute error of 1.8%. More significant deviations (4-5% or 143-148 cm-1) are seen for the NH2 symmetric stretch located on the cytosine moiety as well as for the asymmetric NH2 stretch located on the guanine moiety. This deviation may be due to the smaller basis set used and, as was shown for amino acids, can be improved by using MP2 level potentials. E. Conclusions To date, no simple and accurate method exists which allows the simulation of quantum nuclear dynamics for large systems using PES from high level electronic structure theory. This paper presents the theory and computational implementation of the AICSP method which combines the CSP method for approximate nuclear quantum dynamics with ab initio PES. The method is based on a separable ansatz for the total nuclear wavefunction which results in a set of timedependent single mode equations. Each mode evolves under an average time-dependent potential obtained from classical trajectories simulations for the system, which allows energy to flow between modes. This classical averaging decouples the equations which are implicitly coupled in TDSCF and, as a result, the CSP method is much cheaper to implement than TDSCF. After


24

Page 25 of 33

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


building the single-mode potentials, they are averaged over time and are used to obtain eigenvalues and eigenstates of each mode. The method is applied to stationary vibrational spectroscopy of several amino acids and of HNO3. The agreement obtained between calculated and experimental frequencies is good. The main conclusions of this study are: 1) The level of agreement with experimental results obtained here for several test systems supports the applicability of the AICSP method for stationary vibrational frequencies for other systems as well. 2) The accuracy obtained for AICSP is found to be similar in quality to that of VSCF. Similarly to VSCF, the CSP method provides a fairly accurate description of anharmonic effects without treating them as small perturbations. Unlike most typical applications of VSCF (and TDSCF), where a pairwise approximation is evoked to simplify the calculation of the effective single-mode potentials, AICSP describes many-body effects, albeit classically. 3) AICSP is expected to be efficient for large systems, when high-frequency hydrogen stretching modes are considered, since the required classical trajectories are short. In addition, since the modes are coupled only implicitly, through the time dependence of the average single-mode potentials, an efficient parallelization of the algorithm is possible. It seems very attractive to try and apply the method for hydrogen stretching vibrations in larger molecules, e.g. peptides. 4) An exciting prospect of application of the method is to problems which are fundamentally time-dependent, including vibrational lineshapes and linewidths as well as time-dependent vibrational spectroscopy. AUTHOR INFORMATION Corresponding Author *Email: [email protected]


25


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 26 of 33

ASSOCIATED CONTENT Supporting Information. Tabulated experimental and calculated frequencies are given in the SI as well as an example for the time-dependence of the effective single-mode potential for the OH stretch of conformer I of alanine. This information is available free of charge via the Internet at http://pubs.acs.org

ACKNOWLEDGMENT BH is supported through an Adams fellowship of the Israel Academy of Sciences and Humanities. REFERENCES 1. Hafner, J., Ab-initio Simulations of Materials using VASP: Density-functional Theory and Beyond. J. Comput. Chem. 2008, 29 (13), 2044-78. 2. Gerber, R. B.; Shemesh, D.; Varner, M. E.; Kalinowski, J.; Hirshberg, B., Ab Initio and Semi-empirical Molecular Dynamics Simulations of Chemical Reactions in Isolated Molecules and in Clusters. Phys. Chem. Chem. Phys. 2014, 16, 9760-9775. 3. Ufimtsev, I. S.; Luehr, N.; Martinez, T. J., Charge Transfer and Polarization in Solvated Proteins from Ab Initio Molecular Dynamics. J. Phys. Chem. Lett. 2011, 2 (14), 1789-1793. 4. Lin, I. C.; Seitsonen, A. P.; Tavernelli, I.; Rothlisberger, U., Structure and Dynamics of Liquid Water from ab Initio Molecular Dynamics-Comparison of BLYP, PBE, and revPBE Density Functionals with and without van der Waals Corrections. J. Chem. Theory Comput. 2012, 8 (10), 3902-10. 5. Hassanali, A. A.; Cuny, J.; Verdolino, V.; Parrinello, M., Aqueous Solutions: State of the Art in Ab Initio Molecular Dynamics. Phil. Trans. R. Soc. A 2014, 372, 20120482. 6. Sun, L. P.; Hase, W. L., Born-Oppenheimer Direct Dynamics Classical Trajectory Simulations. In Reviews in Computational Chemistry, Lipkowitz, K. B.; Larter, R.; Cundari, T. R.; Boyd, D. B., Eds. 2003; Vol. 19, pp 79-146. 7. Car, R.; Parrinello, M., Unified Approach for Molecular Dynamics and Densityfunctional Theory. Phys. Rev. Lett. 1985, 55, 2471-2474. 8. Kosloff, R., Time-dependent Quantum-mechanical Methods for Molecular Dynamics. J. Phys. Chem. 1988, 92 (8), 2087-2100. 9. Tal-Ezer, H.; Kosloff, R.; Schaefer, I., New, Highly Accurate Propagator for the Linear and Nonlinear Schrödinger Equation. J. Sci. Comput. 2012, 53 (1), 211-221.


26

Page 27 of 33

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


10. Leforestier, C.; Bisseling, R. H.; Cerjan, C.; Feit, M. D.; Friesner, R.; Guldberg, A.; Hammerich, A.; Jolicard, G.; Karrlein, W.; Meyer, H.-D.; Lipkin, N.; Roncero, O.; Kosloff, R., A Comparison of Different Propagation Schemes for the Time-dependent Schrödinger Equation. J. Comput. Phys. 1991, 94 (1), 59 - 80. 11. Kosloff, D.; Kosloff, R., A Fourier Method Solution for the Time-dependent Schrödinger Equation as a Tool in Molecular Dynamics. J. Comput. Phys. 1983, 52 (1), 35 - 53. 12. Tal-Ezer, H.; Kosloff, R., An Accurate and Efficient Scheme for Propagating the Timedependent Schrödinger Equation. J. Chem. Phys. 1984, 81 (9), 3967. 13. Peskin, U.; Moiseyev, N., The Solution of the Time-dependent Schrödinger Equation by the (t,t’) Method: Theory, Computational Algorithm and Applications. J. Chem. Phys. 1993, 99 (6), 4590. 14. Makri, N., Time-dependent Quantum Methods for Large Systems. Annu. Rev. Phys. Chem. 1999, 50 (1), 167-191. 15. Gerber, R. B.; Ratner, M. A., Mean-field Models for Molecular States and Dynamics: New Developments. J. Phys. Chem. 1988, 92 (11), 3252-3260. 16. Gerber, R. B., Buch, V. , Ranter, M. A., Time-dependent Self-consistent Field Approximation for Intramolecular Energy Transfer. I. Formulation and Application to Dissociation of van der Waals Molecules. J. Chem. Phys. 1982, 77, 3022. 17. Heller, E. J., Time-dependent Variational Approach to Semiclassical Dynamics. J. Chem. Phys. 1976, 64 (1), 63. 18. Alimi, R.; Gerber, R. B., Solvation Effects on Chemical-reaction Dynamics in Clusters Photodissociation of HI in (Xe)NHI. Phys. Rev. Lett. 1990, 64 (12), 1453-1456. 19. García-Vela, A.; Gerber, R. B.; Imre, D. G., Mixed Quantum Wave Packet/Classical Trajectory Treatment of the Photodissociation Process ArHCl→Ar+H+Cl. J. Chem. Phys. 1992, 97 (10), 7242. 20. Gerber, R. B.; Alimi, R., Quantum Effects in Molecular Reaction Dynamics in Solids: Photodissociation of HI in Solid Xe. Chem. Phys. Lett. 1990, 173 (4), 393 - 396. 21. García-Vela, A.; Gerber, R. B., Three-dimensional Quantum Wave Packet Study of the Ar–HCl Photodissociation: A Comparison Between Time-dependent Self-consistent-field and Exact Treatments. J. Chem. Phys. 1995, 103 (9), 3463. 22. Alimi, R.; Gerber, R. B.; Hammerich, a. D.; Kosloff, R.; Ratner, M. a., Validity of Timedependent Self-consistent-field (TDSCF) Approximations for Unimolecular Dynamics: A Test for Photodissociation of the Xe–HI Cluster. J. Chem. Phys. 1990, 93, 6484. 23. Bisseling, R. H.; Kosloff, R.; Gerber, R. B.; Ratner, M. A.; Gibson, L.; Cerjan, C., Exact Time-dependent Quantum Mechanical Dissociation Dynamics of I2He: Comparison of Exact Time-dependent Quantum Calculation with the Quantum Time-dependent Self-consistent field (TDSCF) Approximation. J. Chem. Phys. 1987, 87, 2760. 24. Buch, V.; Gerber, R. B.; Ratner, M. A., Dissociation Dynamics of Ar3 in the Timedependent Self-consistent Field (TDSCF) Approximation. Chem. Phys. Lett. 1983, 101 (1), 44 48. 25. Gerber, R. B.; Alimi, R., Quantum Molecular Dynamics by a Perturbation-Corrected Time-dependent Self-consitent-field Method. Chem. Phys. Lett. 1991, 184 (1), 69 - 75. 26. McCoy, A. B., Transition State Dynamics of Chemical Reactions in Clusters: A Sixdimensional Study of Ar(ClHCl). J. Chem. Phys. 1995, 103 (3), 986. 27. Jungwirth, P.; Gerber, R. B., Quantum Dynamics of Large Polyatomic Systems using a Classically Based Separable Potential Method. J. Chem. Phys. 1995, 102 (15), 6046-6056.


27


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 28 of 33

28. Jungwirth, P.; Fredj, E.; Gerber, R. B., Ultrafast Quantum Dynamics and Resonance Raman Spectroscopy of Photoexcited I2(B) in Large Argon and Xenon Clusters. J. Chem. Phys. 1996, 104, 9332. 29. Jungwirth, P.; Gerber, R. B., Quantum Dynamics of Many-atom Systems by the Classically Based Separable Potential (CSP) Method: Calculations for I−(Ar)12 in Full Dimensionality. J. Chem. Phys. 1995, 102, 8855. 30. Žďánská, P.; Slavı́ček, P.; Jungwirth, P., HCl Photodissociation on Argon Clusters: Effects of Sequential Solvation and Librational Preexcitation. J. Chem. Phys. 2000, 112 (24), 10761. 31. Knospe, O.; Jungwirth, P., Electron Photodetachment in C60−: Quantum Molecular Dynamics with a Non-empirical, `On-the-fly' Calculated Potential. Chem. Phys. Lett. 2000, 317 (6), 529 - 534. 32. Jungwirth, P.; Schmidt, B., Quantum Dynamics Following Electron Photodetachment in the I-Ar2 complex. How Good are the New Separable and Non-separable Simulation Methods? Chem. Phys. Lett. 1997, 275 (3), 127 - 136. 33. Jungwirth, P.; Gerber, R. B., Quantum Molecular Dynamics of Ultrafast Processes in Large Polyatomic Systems. Chem. Rev. 1999, 99 (6), 1583-1606. 34. Makri, N.; Miller, W. H., Time-dependent Self-consistent Field (TDSCF) Approximation for a Reaction Coordinate Coupled to a Harmonic Bath: Single and Multiple Configuration Treatments. J. Chem. Phys. 1987, 87 (1987), 5781. 35. Manthe, U.; Meyer, H. D.; Cederbaum, L. S., Wave-packet Dynamics within the Multiconfiguration Hartree Framework: General Aspects and Application to NOCl. J. Chem. Phys. 1992, 97 (5), 3199. 36. Meyer, H. D.; Manthe, U.; Cederbaum, L. S., The Multi-configurational Time-dependent Hartree Approach. Chem. Phys. Lett. 1990, 165 (1), 73 - 78. 37. Beck, M. H.; Jackle, A.; Worth, G. A.; Meyer, H. D., The Multiconfiguration Timedependent Hartree (MCTDH) Method: A Highly Efficient Algorithm for Propagating Wavepackets. Phys. Rep. 2000, 324 (1), 1 - 105. 38. Meyer, H.-D., Studying Molecular Quantum Dynamics with the Multiconfiguration Time-dependent Hartree Method. WIREs Comput. Mol. Sci. 2012, 2 (2), 351-374. 39. Wang, H.; Thoss, M., Multilayer Formulation of the Multiconfiguration Time-dependent Hartree Theory. J. Chem. Phys. 2003, 119 (3), 1289. 40. Vendrell, O.; Meyer, H. D., Multilayer Multiconfiguration Time-dependent Hartree Method: Implementation and Applications to a Henon-Heiles Hamiltonian and to Pyrazine. J. Chem. Phys. 2011, 134 (4), 044135. 41. Manthe, U., A Multilayer Multiconfigurational Time-dependent Hartree Approach for Quantum Dynamics on General Potential Energy Surfaces. J. Chem. Phys. 2008, 128 (16), 164116. 42. Richter, F.; Hochlaf, M.; Rosmus, P.; Gatti, F.; Meyer, H. D., A Study of the Modeselective Trans--cis Isomerization in HONO using Ab initio Methodology. J. Chem. Phys. 2004, 120 (3), 1306-17. 43. Viel, A.; Eisfeld, W.; Neumann, S.; Domcke, W.; Manthe, U., Photoionization-induced Dynamics of Ammonia: Ab Initio Potential Energy Surfaces and Time-dependent Wave Packet Calculations for the Ammonia Cation. J. Chem. Phys. 2006, 124 (21), 214306.


28

Page 29 of 33

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


44. Richings, G. W.; Polyak, I.; Spinlove, K. E.; Worth, G. A.; Burghardt, I.; Lasorne, B., Quantum Dynamics Simulations using Gaussian Wavepackets: the vMCG Method. Int. Rev. Phys. Chem. 2015, 34 (2), 269-308. 45. Mendive-Tapia, D.; Lasorne, B.; Worth, G. A.; Robb, M. A.; Bearpark, M. J., Towards Converging Non-adiabatic Direct Dynamics Calculations Using Frozen-width Variational Gaussian Product Basis Functions. J. Chem. Phys. 2012, 137 (22), 22A548. 46. Polyak, I.; Allan, C. S.; Worth, G. A., A Complete Description of Tunnelling using Direct Quantum Dynamics Simulation: Salicylaldimine Proton Transfer. J. Chem. Phys. 2015, 143 (8), 084121. 47. Habershon, S.; Manolopoulos, D. E.; Markland, T. E.; Miller, T. F., 3rd, Ring-polymer Molecular Dynamics: Quantum Effects in Chemical Dynamics from Classical Trajectories in an Extended Phase Space. Annu. Rev. Phys. Chem. 2013, 64, 387-413. 48. Cao, J.; Voth, G. A., The Formulation of Quantum Statistical Mechanics Based on the Feynman Path Centroid Density. I. Equilibrium Properties. J. Chem. Phys. 1994, 100 (7), 50935105. 49. Cao, J.; Voth, G. A., The Formulation of Quantum Statistical Mechanics Based on the Feynman Path Centroid Density. II. Dynamical Properties. J. Chem. Phys. 1994, 100 (7), 51065117. 50. Marx, D.; Parrinello, M., Ab Initio Path-integral Molecular Dynamics. Z. Phys. B. 1994, 95, 143-144. 51. Marx, D.; Parrinello, M., Ab Initio Path Integral Molecular Dynamics: Basic Ideas. J. Chem. Phys. 1996, 104 (11), 4077. 52. Tuckerman, M. E.; Marx, D.; Klein, M. L.; Parrinello, M., Efficient and General Aalgorithms for Path Integral Car–Parrinello Molecular Dynamics. J. Chem. Phys. 1996, 104 (14), 5579. 53. Kaczmarek, A.; Shiga, M.; Marx, D., Quantum Effects on Vibrational and Electronic Spectra of Hydrazine Studied by “On-the-Fly” ab Initio Ring Polymer Molecular Dynamics. J. Phys. Chem. A 2009, 113 (10), 1985-1994. 54. Marsalek, O.; Markland, T. E., Ab initio Molecular Dynamics with Nuclear Quantum Effects at Classical Cost: Ring Polymer Contraction for Density Functional Theory. J. Chem. Phys. 2016, 144 (5), 054112. 55. Kapil, V.; VandeVondele, J.; Ceriotti, M., Accurate Molecular Dynamics and Nuclear Quantum Effects at Low Cost by Multiple Steps in Real and Imaginary Time: Using Density Functional Theory to Accelerate Wavefunction Methods. J. Chem. Phys. 2016, 144 (5), 054111. 56. John, C.; Spura, T.; Habershon, S.; Kuhne, T. D., Quantum Ring-Polymer Contraction Method: Including Nuclear Quantum Effects at No Additional Computational Cost in Comparison to Ab Initio Molecular Dynamics. Phys. Rev. E 2016, 93, 043305. 57. Craig, I. R.; Manolopoulos, D. E., Quantum Statistics and Classical Mechanics: Real Time Correlation Functions from Ring Polymer Molecular Dynamics. J. Chem. Phys. 2004, 121 (8), 3368-73. 58. Witt, A.; Ivanov, S. D.; Shiga, M.; Forbert, H.; Marx, D., On the Applicability of Centroid and Ring Polymer Path Integral Molecular Dynamics for Vibrational Spectroscopy. J. Chem. Phys. 2009, 130 (19), 194510. 59. Habershon, S.; Fanourgakis, G. S.; Manolopoulos, D. E., Comparison of Path Integral Molecular Dynamics Methods for the Infrared Absorption Spectrum of Liquid Water. J. Chem. Phys. 2008, 129 (7), 074501.


29


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 30 of 33

60. Miller, W. H., The Semiclassical Initial Value Representation: A Potentially Practical Way for Adding Quantum Effects to Classical Molecular Dynamics Simulations. J. Phys. Chem. A. 2001, 105 (13), 2942--2955. 61. Makri, N.; Miller, W. H., Coherent State Semiclassical Initial Value Representation for the Boltzmann Operator in Thermal Correlation Functions. J. Chem. Phys. 2002, 116 (21), 92079212. 62. Heller, E. J., Time‐dependent Approach to Semiclassical Dynamics. J. Chem. Phys. 1975, 62 (4), 1544-1555. 63. Heller, E. J., Frozen Gaussians: A Very Simple Semiclassical Approximation. J. Chem. Phys. 1981, 75 (6), 2923-2931. 64. Tatchen, J.; Pollak, E., Semiclassical On-the-fly Computation of the S(0)-->S(1) Absorption Spectrum of Formaldehyde. J. Chem. Phys. 2009, 130 (4), 041103. 65. Ceotto, M.; Atahan, S.; Shim, S.; Tantardini, G. F.; Aspuru-Guzik, A., First-principles Semiclassical Initial Value Representation Molecular Dynamics. Phys. Chem. Chem. Phys. 2009, 11 (20), 3861-7. 66. Wehrle, M.; Sulc, M.; Vanicek, J., On-the-fly Ab Initio Semiclassical Dynamics: Identifying Degrees of Freedom Essential for Emission Spectra of Oligothiophenes. J. Chem. Phys. 2014, 140 (24), 244114. 67. Møller, C.; Plesset, M. S., Note on an Approximation Treatment for Many-electron Systems. Phys. Rev. 1934, 46, 618-622. 68. Chai, J. D.; Head-Gordon, M., Long-range Corrected Hybrid Density Functionals with Damped Atom-atom Dispersion Corrections. Phys. Chem. Chem. Phys. 2008, 10 (44), 6615-20. 69. Dirac, P. A. M., Note on Exchange Phenomena in the Thomas Atom. Math. Proc. Cambridge 1930, 26 (03), 376. 70. Lowdin, P. O.; P.K., M., Some Comments on the Time-dependent Variation Principle. Chem. Phys. Lett. 1972, 14 (1), 1 - 7. 71. Chaban, G. M.; Jung, J. O.; Gerber, R. B., Ab Initio Calculation of Anharmonic Vibrational States of Polyatomic Systems: Electronic Structure Combined with Vibrational Selfconsistent Field. J. Chem. Phys. 1999, 111 (5), 1823. 72. Cheng, X.; Steele, R. P., Efficient Anharmonic Vibrational Spectroscopy for Large Molecules using Local-mode Coordinates. J. Chem. Phys. 2014, 141 (10), 104105. 73. Suwan, I.; Gerber, R. B., VSCF in Internal Coordinates and the Calculation of Anharmonic Torsional Mode Transitions. Chem. Phys. 2010, 373 (3), 267-273. 74. Shao, Y.; Gan, Z.; Epifanovsky, E.; Gilbert, A. T. B.; Wormit, M.; Kussmann, J.; Lange, A. W.; Behn, A.; Deng, J.; Feng, X.; Ghosh, D.; Goldey, M.; Horn, P. R.; Jacobson, L. D.; Kaliman, I.; Khaliullin, R. Z.; Ku, T.; Landau, A.; Liu, J.; Proynov, E. I.; Rhee, Y. M.; Richard, R. M.; Rohrdanz, M. A.; Steele, R. P.; Sundstrom, E. J.; Woodcock III, L. H.; Zimmerman, P. M.; Zuev, D.; Albrecht, B.; Alguire, E.; Austin, B.; Beran, G. J. O.; Bernard, Y. A.; Berquist, E.; Brandhorst, K.; Bravaya, K. B.; Brown, S. T.; Casanova, D.; Chang, C.-M.; Chen, Y.; Chien, S. H.; Closser, K. D.; Crittenden, D. L.; Diedenhofen, M.; Distasio, R. A.; Do, H.; Dutoi, A. D.; Edgar, R. G.; Fatehi, S.; Fusti-Molnar, L.; Ghysels, A.; Golubeva-Zadorozhnaya, A.; Gomes, J.; Hanson-Heine, M. W. D.; Harbach, P. H. P.; Hauser, A. W.; Hohenstein, E. G.; Holden, Z. C.; Jagau, T.-C.; Ji, H.; Kaduk, B.; Khistyaev, K.; Kim, J.; Kim, J.; King, R. A.; Klunzinger, P.; Kosenkov, D.; Kowalczyk, T.; Krauter, C. M.; Lao, K. U.; Ele, A.; Laurent, D.; Lawler, K. V.; Levchenko, S. V.; Lin, C. Y.; Liu, F.; Livshits, E.; Lochan, R. C.; Luenser, A.; Manohar, P.; Manzer, S. F.; Mao, S.-P.; Mardirossian, N.; Marenich, A. V.; Maurer, S. A.; Mayhall, N. J.;


30

Page 31 of 33

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


Neuscamman, E.; Oana, C. M.; Olivares-Amaya, R.; O 'neill, D. P.; Parkhill, J. A.; Perrine, T. M.; Peverati, R.; Prociuk, A.; Rehn, D. R.; Rosta, E.; Russ, N. J.; Sharada, S. M.; Sharma, S.; Small, D. W.; Sodt, A.; Stein, T.; Stück, D.; Su, Y.-C.; Thom, A. J. W.; Tsuchimochi, T.; Vanovschi, V.; Vogt, L.; Vydrov, O.; Wang, T.; Watson, M. A.; Wenzel, J.; White, A.; Williams, C. F.; Yang, J.; Yeganeh, S.; Yost, S. R.; You, Z.-Q.; Zhang, I. Y.; Zhang, X.; Zhao, Y.; Brooks, B. R.; Chan, G. K. L.; Chipman, D. M.; Cramer, C. J.; Goddard III, W. A.; Gordon, M. S.; Hehre, W. J.; Klamt, A.; Schaefer III , H. F.; Schmidt, M. W.; Sherrill, C. D.; Truhlar, D. G.; Warshel, A.; Xu, X.; Aspuru-Guzik, A.; Baer, R.; Bell, A. T.; Besley, N. A.; Chai, J.-D.; Dreuw, A.; Dunietz, B. D.; Furlani, T. R.; Gwaltney, S. R.; Hsu, C.-P.; Jung, Y.; Kong, J.; Lambrecht, D. S.; Liang, W.; Ochsenfeld, C.; Rassolov, V. A.; Slipchenko, L. V.; Subotnik, J. E.; Voorhis, T. V.; Herbert, J. M.; Krylov, A. I.; Gill, P. M. W.; Head-Gordon, M., Advances in Molecular Quantum Chemistry Contained in the Q-Chem 4 Program Package. Mol. Phys. 2015, 113, 184-215. 75. Dunning Jr, T. H., Gaussian Basis Sets for Use in Correlated Molecular Calculations. I. The Atoms Boron through Neon and Hydrogen. J. Chem. Phys. 1989, 90, 1007. 76. Herbert, J. M.; Head-Gordon, M., Accelerated, Energy-conserving Born-Oppenheimer Molecular Dynamics via Fock Matrix Extrapolation. Phys. Chem. Chem. Phys. 2005, 7 (18), 3269-75. 77. Barth, A.; Zscherp, C., What Vibrations Tell About Proteins. Q. Rev. Biophys. 2002, 35 (4), 369-430. 78. Barth, A., The Infrared Absorption of Amino Acid Side Chains. Prog. Biophys. Mol. Bio. 2000, 74 (3–5), 141 - 173. 79. Stepanian, S. G.; Reva, I. D.; Radchenko, E. D.; Adamowicz, L., Conformers of Nonionized Proline . Matrix-Isolation Infrared and Post-Hartree - Fock ab Initio Study. J. Phys. Chem. A 2001, 105, 10664-10672. 80. Reva, I. D.; Stepanian, S. G.; Plokhotnichenko, A. M.; Radchenko, E. D.; Sheina, G. G.; Blagoi, Y. P., Infrared Matrix Isolation Studies of Amino Acids. Molecular Structure of Proline. J. Mol. Struct. 1994, 318, 1 - 13. 81. Brauer, B.; Chaban, G. M.; Gerber, R. B., Spectroscopically-tested, Improved, Semiempirical Potentials for Biological Molecules: Calculations for Glycine, Alanine and Proline. Phys. Chem. Chem. Phys. 2004, 6, 2543-2556. 82. Stepanian, S.; Reva, I., Conformational Behavior of Alfa-alanine. Matrix-isolation Infrared and Theoretical DFT and Ab Initio Study. J. Phys. Chem. A 1998, 102, 4623-4629. 83. Stepanian, S. G.; Reva, I. D.; Radchenko, E. D.; Rosado, M. T. S.; Duarte, M. L. T. S.; Fausto, R.; Adamowicz, L., Matrix-Isolation Infrared and Theoretical Studies of the Glycine Conformers. J. Phys. Chem. A 1998, 102, 1041-1054. 84. Kristofer, R. L.; Nathan, P. W.; Keetra, S. P.; and James, A. P., Integrated Intensities of O−H Stretching Bands: Fundamentals and Overtones in Vapor-Phase Alcohols and Acids. J. Phys. Chem. A 2001, 105 (14), 3481-3486. 85. Donaldson, D. J.; Orlando, J. J.; Amann, S.; Tyndall, G. S.; Proos, R. J.; Henry, B. R.; Vaida, V., Absolute Intensities of Nitric Acid Overtones. J. Phys. Chem. A 1998, 102 (27), 51715174. 86. Fleming, P. R.; Li, M.; Rizzo, T. R., Infrared Spectroscopy of Vibrationally Excited HONO2: Shedding Light on the Dark States of Intramolecular Vibrational Energy Redistribution. J. Chem. Phys. 1991, 94 (4), 2425.


31


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 32 of 33

87. Sinha, A.; Vander Wal, R. L.; Crim, F. F., The Vibrationally Mediated Photodissociation Dynamics of Nitric Acid. J. Chem. Phys. 1989, 91 (5), 2929. 88. Miller, Y.; Chaban, G. M.; Gerber, R. B., Theoretical Study of Anharmonic Vibrational Spectra of HNO3, HNO3–H2O, HNO4: Fundamental, Overtone and Combination Excitations. Chem. Phys. 2005, 313 (1–3), 213 - 224. 89. Jung, J. O.; Gerber, R. B., Vibrational Wave Functions and Spectroscopy of (H2O)n, n=2,3,4,5: Vibrational Self-consistent Dield with Correlation Corrections. J. Chem. Phys. 1996, 105 (23), 10332. 90. Abo-Riziq, A.; Grace, L.; Nir, E.; Kabelac, M.; Hobza, P.; de Vries, M. S., Photochemical Selectivity in Guanine-cytosine Base-pair Structures. Proc. Natl. Acad. Sci. USA 2005, 102 (1), 20-3. 91. Brauer, B.; Gerber, R. B.; Kabeláč, M.; Hobza, P.; Bakker, J. M.; Abo Riziq, A. G.; de Vries, M. S., Vibrational Spectroscopy of the G· · ·C Base Pair: Experiment, Harmonic and Anharmonic Calculations, and the Nature of the Anharmonic Couplings. J. Phys. Chem. A 2005, 109 (31), 6974-6984.


32

Page 33 of 33

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


For Table of Contents use only Approximate Quantum Dynamics using Ab Initio Classical Separable Potentials: Spectroscopic Applications Barak Hirshberg, Lior Sagiv and R. Benny Gerber


33

Approximate Quantum Dynamics using Ab Initio ... - ACS Publications

Recommend Documents