VCI, and Ab Initio Molecular - ACS Publications

Subscriber access provided by UNIV OF NEW ENGLAND ARMIDALE

Clusters, Radicals, and Ions; Environmental Chemistry 2

2

IR Spectra of (HCOOH) and (DCOOH): Experiment, VSCF/ VCI and ab Initio Molecular Dynamics Calculations Using Full Dimensional Potential and Dipole Moment Surfaces Chen Qu, and Joel M Bowman J. Phys. Chem. Lett., Just Accepted Manuscript • DOI: 10.1021/acs.jpclett.8b00447 • Publication Date (Web): 30 Apr 2018 Downloaded from http://pubs.acs.org on April 30, 2018

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 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 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.

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 22 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

The Journal of Physical Chemistry Letters

IR Spectra of (HCOOH)2 and (DCOOH)2: Experiment, VSCF/VCI and ab Initio Molecular Dynamics Calculations Using Full Dimensional Potential and Dipole Moment Surfaces Chen Qu and Joel M. Bowman∗ Department of Chemistry and Cherry L. Emerson Center for Scientific Computation, Emory University, Atlanta, Georgia 30322, U.S.A. E-mail: [email protected]

1

ACS Paragon Plus Environment

The Journal of Physical Chemistry Letters 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

Abstract We report quantum VSCF/VCI and ab initio molecular dynamics (AIMD) calculations of the IR spectra of (HCOOH)2 and (DCOOH)2 , using full-dimensional, ab initio potential energy and dipole moment surfaces (PES and DMS). These surfaces are fits, using permutationally invariant polynomials, to 13,475 ab initio CCSD(T)-F12a electronic energies and MP2 dipole moments. Here, “AIMD” means using these ab initio potential and dipole moment surfaces in the MD calculations. The VSCF/VCI calculations use all (24) normal modes for coupling, with a 4-mode representation of the potential. The quantum spectra align well with jet-cooled and room-temperature experimental spectra over the spectral range 600–3,600 cm−1 . Analyses of the complex O-H and C-H stretch bands are made based on the mixing of the VSCF/VCI basis functions. The comparisons of the AIMD IR spectra with both experimental and VSCF/VCI ones provide tests of the accuracy of the AIMD approach. These indicate good accuracy for simple bands but not for the complex O-H stretch band, which is upshifted from experimental and VSCF/VCI bands by roughly 300 cm−1 . In addition to testing the AIMD approach, the PES, DMS and VSCF/VCI calculations for formic acid dimer provide opportunities for testing other methods to represent high-dimensional data and other methods that perform post-harmonic vibrational calculations.

Caclulated IR Spectrum of Formic Acid Dimer

2


Page 2 of 22

Page 3 of 22 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


The formic acid dimer (FAD) (and isotopologues) is one of the mostly intensively studied hydrogen-bonded dimers. It is the smallest example of a carboxylic acid dimer; nevertheless, with 10 atoms, FAD presents a formidable challenge to theory. Recalling that there are two equivalent double hydrogen-bond minima, already signals possibly complex vibrational dynamics. Indeed, this H-bonding motif gives rise to a small tunneling splitting, and possibly contributes, at least indirectly, to some complex bands in the IR (and Raman) spectrum. 1–12 Various aspects of the dynamics of the double-proton transfer and the complex IR (and Raman) spectra have been studied theoretically. 4,10,12–28 Recently, major progress has been made in calculating the tunneling splitting in full dimensionality, 28 using our recent fulldimensional ab initio potential energy surface (PES). 27 The calculated splitting is in good agreement with experimental value of roughly 0.015 cm−1 . 5,11 (For references to the extensive previous literature on calculating the tunneling splitting in reduced dimensionality, see ref. 28.) The present paper is concerned with the IR spectrum. As noted already, the experimental IR (and Raman) spectrum has been reported numerous times, notably in the C-O, C-H and O-H stretch regions of the spectrum. 1–4,9,10,19 Theoretical work on the vibrational spectrum has also been extensive, beginning in 1987 13 when ab initio double-harmonic calculations were reported. Shortly after, numerous reduced-dimensionality calculations were reported. 14,16,18 Recent theoretical work 12,29,30 focused on O-H stretch region of the IR spectrum, owing to the extreme complexity of the experimental band. Florio et al. 12 developed full and reduced dimensionality cubic force fields expanded about a minimum, based on B3LYP calculations. These were used, together with linear dipole moments, to simulate this complex band. With some empirical adjustments to the dipole moment and some force constants, theory could be brought into reasonable accord with experiment. However, the authors noted a an earlier 3 degree-of-freedom (dof) study of the double proton transfer, 18 which reported that excitation of the O-H stretch would result in rapid double-proton transfer. If this were correct, this would lead to a very broad band. We return to this possibility

3



later in the paper. More recently, Matanović and Doˇslić 29 reported a local 2d potential energy surface, based on B3LYP calculations at the minimum and used both second-order perturbative treatment and 2- and 3-mode vibrational configuration interaction calculations. The conclusions were that coupling to low-frequency dimer and rocking modes and Fermi resonances contributed to the width of the band, but had only a minor effect on the frequency of the band maximum. Pitsevich et al. 30 performed VPT2 calculations of the fundamentals of FAD and developed extended 1d and 2d potential energy surfaces using B3LYP/cc-pVTZ calculations to study anharmonic effects in the C-H and O-H stretches. However, the IR spectrum was not simulated in these calculations. All of these calculations are limited in both the dimensionality of the potential, the degree of mode coupling and also the restriction to a single minimum, with the exception of the work by Vener et el. 18 As noted above, we recently reported a full-dimensional ab initio PES, 27 which was used successfully to obtain a small tunneling splitting in excellent agreement with experiment. More recently, we reported a full-dimensional ab initio dipole moment surface (DMS), which was used with the PES in 15 dof VSCF/VCI calculations of the IR spectrum of FAD. 31 Here, we use these potential and dipole moment surfaces in full-dof (24 modes) VSCF/VCI calculations of the IR spectrum of FAD and d-FAD. In addition, we report and assess the accuracy of “ab initio molecular dynamics (AIMD)” calculations of these spectra. We are motivated to perform these calculations for general and specific reasons. The specific motivation is that molecular dynamics (MD) simulations of the IR spectra of (HCOOH)2 and (DCOOH)2 , using model, full-dimensional potential and dipole moment surfaces, were recently reported. 10 In those calculations, the double-proton barrier height and some force constant were adjusted to bring the MD IR spectrum in the complex O-H stretch region into agreement with the room-temperature experimental spectra of FAD and d-FAD. The general motivation is that MD, and more specifically the “AIMD” approach, is widely used to obtain the IR spectra of larger molecules, clusters, etc., owing to difficulties in accurate quantum calculations for such systems.

4


Page 4 of 22

Page 5 of 22 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


Details of the VSCF/VCI calculations are given below; however, we note that these are done with the code MULTIMODE, 32 which uses the exact Watson Hamiltonian together with the efficient n-mode representation of the potential (and dipole moment) surface. 33 These calculations used one minimum as the reference configuration and so they do not describe proton transfer, even though the PES does. (We discuss the possible influence of the double proton transfer, which we conclude is negligible for the present comparison with experiments, later.) The final H-matrix is roughly of total order 80,000, and is 4-fold block diagonal, by making use of the C2h point group of the global minimum. The zero-point energy from this calculation is 15,340 cm−1 , which is in excellent agreement with the value of 15,337±7 cm−1 , reported previously, 27 using the diffusion Monte Carlo method. Thus, we believe the spectrum obtained using the full 24 modes should be reasonably well-converged, at least for the purpose of graphical comparison to experimental spectra. Two sets of MD simulations were done. In one zero-point energy is given initially to all the normal modes. This follows the approach suggested and tested on trans-HONO. 34 The authors of that paper termed this a semi-classically prepared MD but since the method used is the same approximate semi-classical one in quasiclassical trajectory calculations (mainly for gas-phase reaction dynamics) we denoted it as QCMD in a recent experiment/theory paper on small protonated water clusters. 35 Since the procedure is an approximate semiclassical one, it is subject to the well-known “zero-point energy leak.” This can lead to artificial broadening of spectral bands and, owing to the high internal energy (15,583 and 14,156 cm−1 , respectively for (HCOOH)2 and (DCOOH)2 ), other unphysical effects, noted below. However, the hope is that it does capture the effects of anharmonicty. The second set of trajectories differs in that the total energy given to the molecule equals 24kB T , so that on average each normal mode has an energy of kB T , i.e., the average classical energy of a harmonic oscillator at temperature T . Large batches of trajectories were run for roughly 12 ps, with initial conditions generated using “normal-mode sampling.” 36 From these the dipole-dipole correlation was obtained and Fourier transformed to obtain the IR spectrum.

5



Page 6 of 22

More details of these calculations are given below and also in the SI. The goal of these MD calculations is to locate the spectral bands, with a particular focus on the complex highfrequency band. Delving into an interpretation of the structure of this band is difficult using the MD approach, and so we do not attempt to do that. The VSCF/VCI and experimental IR spectra of (HCOOH)2 and (DCOOH)2 are shown in Figure 1 and 2, respectively. For reference, the double-harmonic stick spectra are also given; these are from the PES and DMS. Slightly reduced dimensionality, 21-mode spectra are also given. These excluded the three lowest-frequency out-of-plane modes. The normal modes are depicted and labeled in Table S1 of Supporting Information (SI). Figs. S1-S3 in the SI gives further comparisons between the 21-mode and 24-mode spectra as well as 17, 15 and 12-mode spectra. The 21 and 24-mode calculated spectra are in good agreement with each other, with the latter somewhat upshifted from the former. The experimental survey spectra are room temperature ones, 10 whereas the spectra in the high-frequency region are jet-cooled ones, 3 reproduced by digitizing the original plot in that reference. The calculated spectra are for total angular momentum zero and with no hot bands. In these plots the calculated stick intensities are Gaussian-broadened with full width at half maximum (FWHM) of 23.5 cm−1 , in order to achieve a similar width of the room-temperature experiment. In the comparison with the jet-cooled spectra the sticks are similarly broadened but with FWHM of 11.8 cm−1 . Overall, the experimental spectra are reproduced with unprecedented accuracy, notably for the complex band in the range 2,600–3,400 cm−1 The assignment of this complex band is clearly problematic. The VSCF/VCI calculations provide the CI coefficients of the basis for each molecular eigenstate, ψL , from the usual equation, ψL =

N X

(L)

ci χ i ,

(1)

i=1

where the basis, {χ}, consists of direct-product, virtual states of the ground-state VSCF Hamiltonian and N is the size of the basis, which is roughly 20,000 for each symmetry block. For the purposes of discussion, we can assume that these are very similar to the familiar 6


Page 7 of 22 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


standard harmonic-oscillator direct product basis. While the vast majority are very small in magnitude, i.e., less than 0.01; for this complex band, all of the bright small sticks (which are subsequently broadened) have CI coefficients that are at least 0.01 in magnitude for the zero-order O-H or C-H fundamental. In the usual approach, we examine the significant CI coefficients for each bright state and thereby make the “assignment”. This was done in a table in our paper reporting the DMS and VSCF/VCI calculations, which however, just considered 15 in-plane modes. 31 In the present 24 (and 21)-mode calculations, such a tabulation is not illuminating, except to note that these coefficients are all very small. Instead, we present a graphical representation of this analysis for the 24- and 21-mode calculations. Specifically, in Figure 3, we plot the sum of the squares of the VCI coefficients of the C-H and O-H stretches in 20 cm−1 -windows, as a function of energy across the band. One can immediately see that the O-H stretch coefficients are widely spread all over the band. Therefore, it’s essentially impossible to assign a certain range to the O-H stretch. On the other hand, though there are also several states with significant coefficient for the C-H stretch, it is less diluted, and most of them cluster between 2,900 and 3,000 cm−1 , so we could assign this range to the C-H stretch, which agrees with the conclusion from isotopic substitution experiment. 2 We also see good agreement between the 21 and 24-mode results. Clearly, the zeroth-order O-H stretch basis function is highly fractionated among the molecular eigenstates. We can get a quantitative measure of this by considering the inverse expansion of the zeroth-order O-H fundamental and C-H fundamental in terms of molecular eigenstates. From the inverse expansion we determine that there are 272 molecular eigenstates between 2,500 and 3,500 cm−1 that account for roughly 90% of the zeroth-order C-H stretch fundamental , and 522 molecular eigenstates that account for 80% of the zeroth-oder O-H stretch fundamental. Also, it is worth noting that by performing the full 24-mode calculations, we discovered the importance of the out-of-plane modes. One may speculate that both C-H and O-H stretches are in-plane modes and there is not much coupling between those two and the out-of-plane modes, and indeed this separation of modes has been done

7



previously by others 25 and also recently by us. 31 However, our 17-in-plane-mode calculation significantly deviates from the experiment (see Fig. S2 in SI), with a very intense transition at about 2,800 cm−1 . Only when the C-H and O-H out-of-plane bending modes are coupled in the 21 and 24-mode calculations we see the good agreement with the experiment. The large degree of coupling in this spectral range need more discussion, but before doing, we present the MD spectra, to assess the accuracy and also to perhaps gather more insight into the coupling. These spectra for FAD and d-FAD and comparisons with experiment 10 are given with Figure 4. As seen, both the classical NVE and QCMD spectra agree well with experiment for the simple bands that can be charactered as relatively pure fundamentals of normal modes. The complexity of the O-H stretch band is reproduced by the MD simulations; this is a major improvement over the double-harmonic stick spectrum. However, the MD bands are substantially upshifted compared to the experiments and the QM calculations. The center of the broad MD bands is close to the harmonic frequency of the bright OHstretch, 3326 cm−1 . This result for the MD spectra is not surprising, since it is well-known that the classical MD can not reproduce the full anharmonicity, especially for the highfrequency HX stretches. When the temperature is low, the trajectory basically just samples the harmonic part of the PES, and thus the MD spectrum is close to the double-harmonic one, which does not capture the complexity of the low-temperature spectrum in the O-H stretching region. This is shown in Figure S4 , which is the IR spectrum calculated using classical NVE trajectories at 50 K, which is approximately the temperature of the jet-cooled experiment. The hoped-for improvement of the QCMD over the standard classical spectrum is largely not realized. There is some improvement overall; however, for the complex band the most noticeable difference from the classical spectrum is the breath of the band, which is in somewhat better agreement with experiment. The classical MD results, especially for the complex band, are significant in several respects. First, they relate directly to an assessment of room-temperature NVE MD calculations of the IR spectra of FAD and d-FAD using a flexible, model potential energy surface

8


Page 8 of 22

Page 9 of 22 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 a fixed-charge dipole moment surface, mentioned above. 10 In that interesting study, the barrier height of the model potential was varied with the goal of improving agreement with experiment in the complex O-H stretch region of the spectrum. The model potential denoted MMPT-MP2 has a double proton barrier of 8.2 kcal/mol. However, the calculated O-H stretch band was up-shifted from experiment by roughly 300 cm−1 , (in good agreement with the present MD result) and so adjustments to this barrier were made, i.e., lowering it, which resulted in better agreement with experiment. A barrier of 5.4 kcal/mol produced the best agreement with experiment for this band. Based on this and other results where the O-H stretch frequency was lowered, these authors concluded that the correct barrier is between 5 and 7 kcal/mol. Clearly, based on the present results we disagree with this conclusion because the classical MD spectrum for this complex band is not accurate. Further, the barrier of the ab initio PES used here 27 is 8.17 kcal/mol, which is slightly below the current benchmark value of 8.32 kcal/mol reported by Marx and co-workers, 26 based on CCSD(T)/aV5Z calculations. This PES used in the quantum VSCF/VCI calculations of the IR spectra as well recent as instanton calculations of the tunneling splitting, 28 produce very good agreement with experiment and so the barrier of ca. 8.2 kcal/mol is indeed accurate. These recent NVT MD calculations of the FAD and d-FAD spectra are very insightful, as the complexity of the O-H stretch band was reproduced in those calculations and the position of the band was shown to depend sensitively on the barrier height to double-proton transfer. This leads us to the second discussion point about how the classical MD spectra provide additional insight to the major issue of interpreting the complex O-H stretch band. One possibility, raised in calculations by Vener, K¨ uhn and Bowman in 1971, 18 is that excitation of the O-H stretch fundamental results in a large double-well splitting of 70 cm−1 , and thus significant delocalization of the wavefunction over the two wells. That conclusion was based on calculations in 3 dof, so a very reduced dimensionality treatment. The present VSCF/VCI calculations do not support this conclusion, as the vibrational basis used does not span the second well. This localization of the wavefunctions is also consistent with the

9



QCMD calculations, where the total energy of FAD is roughly 4 times the barrier of the double-proton transfer and yet only a small number (about 10%) of the hundred trajectories visited the second well; this occurred because of “ZPE leak”. These few trajectories were not used in the calculation of the IR spectra and so the complexity seen in the QCMD spectra is not due to trajectory motion that visits both wells. Finally, we note the 7mode wavepacket calculations by Luckhaus, 25 who determined the splitting for vibrationally excited states of FAD and reported values much less than 1 cm−1 , which clearly indicates localization of the wavefunctions. So, at the level of resolution of both theory and experiment shown here, neglecting double-proton transfer is justified. Nevertheless, the analysis of the VCI coefficients and examination of QCMD trajectories do indicate that there is significant motion along the direction of the double proton transfer. This intermolecular motion causes further broadening of the band that is strongly fractionated by intramolecular couplings. Taken together, both sets of modes are responsible for the complexity of this band. Turning to DCOOH2 , We note a similar degree fractionization of the O-H stretch, in agreement with experiment. This is not surprising as the double H-bonded O-H stretch is also present in this isotopolog. There is a new feature in this spectrum at around 2,200 cm−1 , which is nominally the C-D stretch. In fact, this band is not a single stick, but a classic Fermi resonance doublet that was reported experimentally and analyzed in the usual way. 6,9 The 24-mode VSCF/VCI do find this doublet due to the near resonance of the C-D stretch and the combination band of the C-O single-bond stretch and C-D in-plane bend, in agreement with experiment. However, the VSCF/VCI calculations show a small intensity for the upper component of this doublet. Interestingly, in the smaller 12 and 15-mode calculation, this feature can be well reproduced (see Fig. S3 in SI). Also, the doublet band is seen in classical spectrum, perhaps fortuitously. We conclude with several comments about possible future work. One is about the fully deuterated FAD, (DCOOD)2 . The Raman spectrum shows that the complexity of the band in the stretching region is substantially narrowed, though not completely gone. 4 The IR

10


Page 10 of 22

Page 11 of 22 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


spectrum should have similar features, and we plan to calculate this spectrum using the present approach to see if this is the case. A second comment relates to theoretical methodology. FAD, with 10 atoms, is currently the largest complex to have full dimensional ab initio potential and dipole moments surface. The data set of energies and dipole moments is of modest size, i.e., around 13000; however, the dimensionality of the space is large. In terms of Morse-transformed internuclear distances, used for the fits, the dimensionality is 45. These data sets are available for fitting by other methods from “machine learning” codes, such as Neural Networks and Gaussian Process Regression, both of which are available in popular codes such as Matlab 37 (both methods) and R 38 (Gaussian Process Regression). The calculation of the complex spectra of FAD is clearly a challenge for theory and fortunately the Watson Hamiltonian, used in the present 24-mode VSCF/VCI calculations, does appear to be suited for the types of motion relevant to these spectra. Other approaches that are extensions of the MD approach to obtain complex spectra, 39,40 such Centroid Molecular Dynamics, 41 Ring Polymer Molecular Dynamics, 42 and the Initial Value Representation 43 could now be usefully tested with the availability of the potential and dipole moment surfaces.

Computational details The VSCF/VCI calculations were performed with the MULTIMODE software, 44 which determines the eigenvalues and eigenfunctions of the Watson Hamiltonian. 45 To make the calculations feasible, a 4-mode representation 33 of the FAD potential was used for the calculation of the eigenvalues, and a 3-mode representation of the dipole moment was used in order to evaluate the dipole matrix elements numerically efficiently. These calculations were similar to the ones we recently reported, 31 except that we coupled more normal modes, up to 24 (all) in the current calculations. The C2h point-group symmetry was exploited to set up a 4-block Hamiltonian matrix, each of which is of order 20,000. Numerous standard convergence tests, including the n-mode representation of the potential, the size of basis

11



functions, and the excitation space in the VCI calculation, were performed and we estimate that the calculated spectra (for a given number of coupled modes) shown are converged to within roughly 10 cm−1 or less. Using 5-mode representation of the potential or using 2-mode representation of dipole produces virtually the same spectrum, so the 4-mode representation of the potential and the 3-mode representation of the dipole are converged in terms of this n-mode expansion. As noted above, two sets of trajectories were propagated, using the present PES and DMS, for roughly 12 ps, which corresponds to 200,000 time steps. In the QCMD set six hundred trajectories were performed. Owing to the high energy of these trajectories, a significant number result in unphysical dimer dissociation or rapid double-proton transfer. These trajectories were discarded and finally 293 QCMD trajectories were used to compute the IR spectrum of (HCOOH)2 and 439 for (DCOOH)2 . Tests indicate that the spectra are well converged with this number of trajectories. In the second set of trajectories the molecule was given a total energy of 24kB T , the average thermal energy at temperature T , for this 24-mode molecule. The total energy is distributed among the normal modes, and then the initial conditions were sampled using the procedures described in ref. 36. In the calculations shown above T = 300 K, and results for T = 50 K are given in the SI. For this set of calculations one hundred trajectories were performed for the given T and again tests indicate that the spectra are well converged with this number. The dipole-dipole correlation function C(t) = µ(t)·µ(0) was obtained for each set of trajectories, then averaged over the set of trajectories, and finally Fourier transformed, without any filtering, to obtain the spectra. The Fourier transformation was done using the “Discrete Fourier Transform” function in Numpy. 46 Tests of the dependence of the final spectrum with respect to the choice of time origin are given in SI. These show virtually no dependence on choice for two cases, one where the origin is the initial time (the correlation function starts at t = 0, without “equilibration”) for the propagation of trajectories and the other where the molecule was equilibrated for 50,000 time steps, or approximately 3 ps, before the dipole-dipole correlation

12


Page 12 of 22

Page 13 of 22 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


function was computed. Finally, in the SI we plot the distribution of kinetic energy vs time for one trajectory and the distribution of the instantaneous temperature in all the trajectories corresponding to T = 300 K. As expected, the average kinetic energy is 12kB T (half of the total energy) and the instantaneous “temperature” obtained from the fluctuations of the kinetic energy is Gaussian with a peak at 300 K and full-width at half-maximum of approximately 200 K. These results are consistent with a microcanonical NVE distribution, as expected.

Acknowledgement We thank Henrik G. Kjaergaard for sending the room temperature experimental spectra and a reviewer for excellent suggestions. We also thank the National Science Foundation CHE-1463552 for financial support.

Supporting Information Available This contains the visualization of the normal modes of FAD, calculated IR spectra of (HCOOH)2 and (DCOOH)2 that are not shown in the main article, and a figure showing the kinetic energy and instantaneous temperature in classical molecular dynamic calculations. This material is available free of charge via the Internet at http://pubs.acs.org/.

References (1) Ito, F.; Nakanaga, T. A Jet-Cooled Infrared Spectrum of the Formic Acid Dimer by Cavity Ring-Down Spectroscopy. Chem. Phys. Lett. 2000, 318, 571–577. (2) Ito, F.; Nakanaga, T. Jet-Cooled Infrared Spectra of the Formic Acid Dimer by Cavity Ring-Down Spectroscopy: Observation of the O-H Stretching Region. Chem. Phys. 2002, 277, 163–169. 13



(3) Georges, R.; Freytes, M.; Hurtmans, D.; Kleiner, I.; Vander Auwera, J.; Herman, M. Jet-Cooled and Room Temperature FTIR Spectra of the Dimer of Formic Acid in the Gas Phase. Chem. Phys. 2004, 305, 187–196. (4) Zielke, P.; Suhm, M. A. Raman Jet Spectroscopy of Formic Acid Dimers: Low Frequency Vibrational Dynamics and Beyond. Phys. Chem. Chem. Phys. 2007, 9, 4528– 4534. ¨ Havenith, M. High-Resolution Infrared Spectroscopy of the Formic Acid (5) Birer, O.; Dimer. Annu. Rev. Phys. Chem. 2009, 60, 263–275. (6) Yoon, Y. H.; Hause, M. L.; Case, A. S.; Crim, F. F. Vibrational Action Spectroscopy of the C-H and C-D Stretches in Partially Deuterated Formic Acid Dimer. J. Chem. Phys. 2008, 128, 084305. (7) Xue, Z.; Suhm, M. A. Probing the Stiffness of the Simplest Double Hydrogen Bond: The Symmetric Hydrogen Bond Modes of Jet-cooled Formic Acid Dimer. J. Chem. Phys. 2009, 131, 054301. (8) Xue, Z.; Suhm, M. A. Adding more weight to a molecular recognition unit: the lowfrequency modes of carboxylic acid dimers. Mol. Phys. 2010, 108, 2279–2288. (9) Nydegger, M. W.; Rock, W.; Cheatum, C. M. 2D IR Spectroscopy of the C-D Stretching Vibration of the Deuterated Formic Acid Dimer. Phys. Chem. Chem. Phys. 2011, 13, 6098–6104. (10) Mackeprang, K.; Xu, Z. H.; Maroun, Z.; Meuwly, M.; Kjaergaard, H. G. Spectroscopy and Dynamics of Double Proton Transfer in Formic Acid Dimer. Phys. Chem. Chem. Phys. 2016, 18, 24654–24662. (11) Zhang, Y.; Li, W.; Luo, W.; Zhu, Y.; Duan, C. High Resolution Jet-cooled Infrared

14


Page 14 of 22

Page 15 of 22 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


Absorption Spectra of (HCOOH)2 , (HCOOD)2 , and HCOOH-HCOOD Complexes in 7.2 µm Region. J. Chem. Phys. 2017, 146, 244306. (12) Florio, G. M.; Zwier, T. S.; Myshakin, E. M.; Jordan, K. D.; Sibert, E. L. Theoretical Modeling of the OH Stretch Infrared Spectrum of Carboxylic Acid Dimers Based on First-Principles Anharmonic Couplings. J. Chem. Phys. 2003, 118, 1735–1746. (13) Chang, Y.-T.; Yamaguchi, Y.; Miller, W. H.; Schaefer, H. F. An Analysis of the Infrared and Raman Spectra of the Formic Acid Dimer (HCOOH)2 . J. Am. Chem. Soc. 1987, 109, 7245–7253. (14) Shida, N.; Barbara, P. F.; Almlöf, J. E. A Reaction Surface Hamiltonian Treatment of the Double Proton Transfer of Formic Acid Dimer. J. Chem. Phys. 1991, 94, 3633. (15) Kim, Y. Direct Dynamics Calculation for the Double Proton Transfer in Formic Acid Dimer. J. Am. Chem. Soc. 1996, 118, 1522–1528. (16) Loerting, T.; Liedl, K. R. Toward Elimination of Discrepancies between Theory and Experiment: Double Proton Transfer in Dimers of Carboxylic Acids. J. Am. Chem. Soc. 1998, 120, 12595–12600. (17) Miura, S.; Tuckerman, M. E.; Klein, M. L. An ab Initio Path Integral Molecular Dynamics Study of Double Proton Transfer in the Formic Acid Dimer. J. Chem. Phys. 1998, 109, 5290–5299. (18) Vener, M. V.; K¨ uhn, O.; Bowman, J. M. Vibrational Spectrum of the Formic Acid Dimer in the OH Stretch Region. A Model 3D Study. Chem. Phys. Lett. 2001, 349, 562–570. (19) Ushiyama, H.; Takatsuka, K. Successive Mechanism of Double-proton Transfer in Formic Acid Dimer: A Classical Study. J. Chem. Phys. 2001, 115, 5903–5912.

15



(20) Markwick, P. R. L.; Doltsinis, N. L.; Marx, D. Targeted Car-Parrinello Molecular Dynamics: Elucidating Double Proton Transfer in Formic Acid Dimer. J. Chem. Phys. 2005, 122, 054112. (21) Fillaux, F. Quantum Entanglement and Nonlocal Proton Transfer Dynamics in Dimers of Formic Acid and Analogues. Chem. Phys. Lett. 2005, 408, 302–306. (22) Luckhaus, D. Concerted Hydrogen Exchange Tunneling in Formic Acid Dimer. J. Phys. Chem. A 2006, 110, 3151–3158. (23) Barnes, G. L.; Sibert, E. L. An Equilibrium Focused Approach to Calculating the Raman Spectrum of the Symmetric OH Stretch in Formic Acid Dimer. J. Mol. Spectrosc. 2008, 249, 78–85. (24) Burisch, C.; Markwick, P. R. L.; Doltsinis, N. L.; Schlitter, J. ’Dynamic Distance’ Reaction Coordinate for Competing Bonds: Applications in Classical and ab Initio Simulations. J. Chem. Theory Comput. 2008, 4, 164–172. (25) Luckhaus, D. Hydrogen Exchange in Formic Acid Dimer: Tunnelling above the Barrier. Phys. Chem. Chem. Phys. 2010, 12, 8357–8361. (26) Ivanov, S. D.; Grant, I. M.; Marx, D. Quantum Free Energy Landscapes from ab Initio Path Integral Metadynamics: Double Proton Transfer in the Formic Acid Dimer is Concerted but not Correlated. J. Chem. Phys. 2015, 143, 124304. (27) Qu, C.; Bowman, J. M. An ab Initio Potential Energy Surface for the Formic Acid Dimer: Zero-point Energy, Selected Anharmonic Fundamental Energies, and Groundstate Tunneling Splitting Calculated in Relaxed 1–4-mode Subspaces. Phys. Chem. Chem. Phys. 2016, 18, 24835. (28) Richardson, J. O. Full- and Reduced-dimensionality Instanton Calculations of the Tun-

16


Page 16 of 22

Page 17 of 22 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


nelling Splitting in the Formic Acid Dimer. Phys. Chem. Chem. Phys. 2017, 19, 966– 970. (29) Matanović, I.; Doˇslić, N. Theoretical Modeling of the Formic Acid Dimer Infrared Spectrum: Shaping the O−H Stretch Band. Chem. Phys. 2007, 338, 121–126. (30) Pitsevich, G. A.; Malevich, A. E.; Kozlovskaya, E. N.; Doroshenko, I. Y.; Sablinskas, V.; Pogorelov, V.; Dovgal, D.; Balevicius, V. Anharmonic Analysis of C-H and O-H Stretching Vibrations of the Formic Acid Dimer. Vib. Spectrosc. 2015, 79, 67–75. (31) Qu, C.; Bowman, J. M. High-dimensional fitting of sparse datasets of CCSD(T) electronic energies and MP2 dipole moments, illustrated for the formic acid dimer and its complex IR spectrum. J. Chem. Phys. 2018, 148, 241713. (32) Huang, X.; Braams, B. J.; Bowman, J. M. Ab initio potential energy and dipole moment surfaces for H5 O2+ . J. Chem. Phys 2005, 122, 044308. (33) Carter, S.; Culik, S. J.; Bowman, J. M. Vibrational Self-consistent Field Method for Many-mode Systems: A New Approach and Application to the Vibrations of CO Adsorbed on Cu(100). J. Chem. Phys. 1997, 107, 10458–10469. (34) Van-Oanh, N.-T.; Falvo, C.; Calvo, F.; Lauvergnat, D.; Basire, M.; Gaigeot, M.P.; Parneix, P. Improving anharmonic infrared spectra using semiclassically prepared molecular dynamics simulations. Phys. Chem. Chem. Phys. 2012, 14, 2381–2390. (35) Esser, T. K.; Knorke, H.; Asmis, K. R.; Schöllkopf, W.; Yu, Q.; Qu, C.; Bowman, J. M.; Kaledin, M. Deconstructing Prominent Bands in the Terahertz Spectra of H7 O3+ and H9 O4+ : Intermolecular Modes in Eigen Clusters. J. Phys. Chem. Letts. 2018, 9, 798– 803. (36) Hase, W. L. In Encyclopeida of Computational Chemistry; Allinger, N. L., Ed.; Wiley: New York, 1998; Vol. 1; pp 399–407. 17



(37) MATLAB, Release R2017b; The MathWorks Inc.: Natick, Massachusetts, 2017. (38) R Core Team, R: A Language and Environment for Statistical Computing. R Foundation for Statistical Computing: Vienna, Austria, https://www.R-project.org (accessed Apr 26, 2018). (39) Rossi, M.; Liu, H.; Paesani, F.; Bowman, J.; Ceriotti, M. Communication: On the consistency of approximate quantum dynamics simulation methods for vibrational spectra in the condensed phase. J. Chem. Phys. 2014, 141, 181101. (40) Ivanov, S. D.; Witt, A.; Shiga, M.; Marx, D. Communications: On artificial frequency shifts in infrared spectra obtained from centroid molecular dynamics: Quantum liquid water. J. Chem. Phys. 2010, 132, 031101. (41) Jang, S.; Voth, G. A. A derivation of centroid molecular dynamics and other approximate time evolution methods for path integral centroid variables. J. Chem. Phys. 1999, 111, 2371–2384. (42) Rossi, M.; Ceriotti, M.; Manolopoulos, D. E. How to remove the spurious resonances from ring polymer molecular dynamics. J. Chem. Phys. 2014, 140, 234116. (43) Liu, J.; Miller, W. H. Using the thermal Gaussian approximation for the Boltzmann operator in semiclassical initial value time correlation functions. J. Chem. Phys. 2006, 125, 224104. (44) Bowman, J. M.; Carter, S.; Huang, X. MULTIMODE: A Code to Calculate Rovibrational Energies of Polyatomic Molecules. Int. Rev. Phys. Chem. 2003, 22, 533. (45) Watson, J. K. G. Simplification of the Molecular Vibration-Rotation Hamiltonian. Mol. Phys. 1968, 15, 479. (46) Oliphant, T. E. Guide to NumPy, 2nd ed.; Continuum Press: Austin, TX, 2015.

18


Page 18 of 22

Page 19 of 22 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


Figure 1: Comparison between the calculated and the experimental IR spectra of (HCOOH)2 : (a) 24-mode calculation vs. room-temperature spectrum and CO, C-H and O-H stretch normal modes; (b) 21-mode calculation vs. jet-cooled spectrum; (c) 24-mode calculation vs. jet-cooled spectrum. 3

19



(a)

24-mode Experiment Double harmonic 1000

1500

2000

2500

3000

3500

(b)

Intensity (arb. unit)

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 22

21-mode Experiment 2850 2900 2950 3000 3050 3100 3150 3200 3250 (c)

24-mode Experiment 2850 2900 2950 3000 3050 3100 3150 3200 3250 Frequency (cm-1)

Figure 2: Comparison between the calculated and the experimental IR spectra of (DCOOH)2 : (a) 24-mode calculation vs. room-temperature spectrum ; (b) 21-mode calculation vs. jetcooled spectrum; 2 (c) 24-mode calculation vs. jet-cooled spectrum.

20


Page 21 of 22

0.4

Coefficient squared

0.35

CH stretch OH stretch

(a) 21-mode

0.3 0.25 0.2 0.15 0.1 0.05 0 2600

2800

3000

3200

3400

0.45 0.4 Coefficient squared

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


CH stretch OH stretch

(b) 24-mode

0.35 0.3 0.25 0.2 0.15 0.1 0.05 0 2600

2800

3000

3200

3400

-1

Frequency (cm )

Figure 3: Square of the VCI coefficient of the C-H and OH stretches as a function of energy across the stretching band, based on (a) 21-mode VSCF/VCI calculation; (b) 24-mode VSCF/VCI calculation.

21



(a)


MD 300K QCMD Experiment

600

1000

1400

1800

2200

2600

(b)

3000

3400

MD 300K QCMD Experiment


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 22

600

1000

1400

1800

2200

2600

3000

3400

-1

Frequency (cm )

Figure 4: Comparison between the spectra calculated by MD simulations and the experiments: (a) classical 300K MD, QCMD, and room-temperature spectrum for (HCOOH)2 ; (b) classical 300K MD, QCMD, and room-temperature spectrum for (DCOOH)2 .

22


VCI, and Ab Initio Molecular - ACS Publications

Recommend Documents