Anharmonic Rovibrational Partition Functions for Fluxional Species at

Subscriber access provided by READING UNIV

Article

Anharmonic Rovibrational Partition Functions for Fluxional Species at High Temperatures via Monte Carlo Phase Space Integrals Ahren W Jasper, Zackery Boone Gruey, Lawrence B Harding, Yuri Georgievskii, Stephen J Klippenstein, and Albert F. Wagner J. Phys. Chem. A, Just Accepted Manuscript • DOI: 10.1021/acs.jpca.7b11722 • Publication Date (Web): 22 Jan 2018 Downloaded from http://pubs.acs.org on January 22, 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 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.

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

J. Phys. Chem. A January 19, 2018

Anharmonic Rovibrational Partition Functions for Fluxional Species at High Temperatures via Monte Carlo Phase Space Integrals Ahren W. Jasper,* Zackery B. Gruey, Lawrence B. Harding, Yuri Georgievskii, Stephen J. Klippenstein, and Albert F. Wagner Chemical Sciences and Engineering Division, Argonne National Laboratory, Argonne, IL 60439

Abstract. Monte Carlo phase space integration (MCPSI) is used to compute full dimensional and fully anharmonic––but classical––rovibrational partition functions for 22 small- and medium-sized molecules and radicals. Several of the species considered here feature multiple minima and low-frequency nonlocal motions, and efficiently sampling these systems is facilitated using curvilinear (stretch, bend, torsion) coordinates. The curvilinear coordinate MCPSI method is demonstrated to be applicable to the treatment of fluxional species with complex rovibrational structures and as many as 21 fully coupled rovibrational degrees of freedom. Trends in the computed anharmonicity corrections are discussed. For many systems, rovibrational anharmonicities at elevated temperatures are shown to vary consistently with the number of degrees of freedom and with temperature once rovibrational coupling and torsional anharmonicity are accounted for. Larger corrections are found for systems with complex vibrational structures, such as systems with multiple large amplitude modes and/or multiple minima.

*E-mail: [email protected]

ACS Paragon Plus Environment

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

2 1. Introduction The neglect of vibrational anharmonicity can be the dominant source of uncertainty in a priori predictions of thermochemistry and kinetics.1,2 Errors in the zero point energy alone can be significant.3,4 In a recent study,5 anharmonic zero point energy corrections as large as 300 cm–1 were calculated using diffusion Monte Carlo6 for 48 small and moderately-sized radicals and molecules. Earlier studies of the highly fluxional ions CH5+ and H5O2+ found that anharmonic zero point energies, again computed using diffusion Monte Carlo, differed from harmonic values by as much as 450 cm–1.7–9 At room temperature, an error of just 150 cm–1 in the threshold or reaction energy leads to an uncertainty of a factor of two in the computed rate coefficient or equilibrium constant. At higher temperatures, a priori accuracies are less sensitive to anharmonicities in the zero point energy, and, instead, rovibrational anharmonicity corrections to state counts and partition functions become significant. These errors are often the major source of uncertainty in a priori calculations at high temperatures. The study in Ref. 5 included a large number of fluxional species, i.e., species with nonlocal motions such as inversions and torsions. Fluxional species were shown to have anharmonic zero point energy corrections approximately twice as large, on average, as those for nonfluxional species. Also demonstrated in Ref. 5, strategies for computing anharmonic zero point energies that rely on a local expansion about a single minimum, such as vibrational perturbation theory 10 (VPT2), were found to be inconsistent as predictive methods for fluxional species. Specifically, errors in the VPT2 calculations were found to be overly sensitive to the choice of basis set for some systems, while for others the error in the VPT2 calculation was larger than the anharmonicity correction itself.5 Here we consider the calculation of anharmonic rovibrational properties (state counts, state densities, and partition functions) at high temperatures and energies, again with an emphasis on the treatment of fluxional species. A variety of methods for computing rovibrational molecular properties have been developed. Exact quantum mechanical approaches include vibrational self consistent field, configuration interaction, and coupled cluster theories,

11 – 21

path integral sampling, 22 – 28 and direct product

diagonalization. 29 – 39 Significant (and typically prohibitive) challenges remain in the


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


3 application of these approaches to systems with more than a few degrees of freedom, nonzero rotation, and high internal energies. Exact methods are often developed for use in spectroscopic studies of low-energy transitions, and they are not frequently applied in the computation of high temperature rovibrational partition functions. More approximate and more broadly applicable approaches include empirical corrections for stretches, bends, and torsions, 40,41 perturbation theory,10,42–48 specialized treatments for torsions, 49–58 and reduced- and one-dimensional anharmonic models.59 While these approaches can be of great practical utility, the validity of the assumptions and empiricism involved in their specific implementations has not been well studied, particularly for fluxional species at high temperatures and energies, with some exceptions.48,56,60 Here we use a semiclassical approach based on classical Monte Carlo phase space integration61–66 (MCPSI). We demonstrate that this method is applicable to fluxional chemical systems at high energies and temperatures and for nonzero overall rotation. Our test set consists of 22 of the 48 species considered in Ref. 5 and includes species with challenging vibrational structures, including torsions, inversions, and strongly coupled modes. The MCPSI approach, while classical, makes no additional assumptions about the nature of nonlocal motions, requires no artificial mode separations, and is not based on an expansion about a preferred minimum. Our implementation67,68 of MCPSI makes use of sampling in curvilinear coordinates, which facilitates the treatment of fluxional species, particularly those with multiple minima. Furthermore, Monte Carlo methods trivially and near ideally exploit the massive numbers of cores that can be made available to calculations. Importantly, the present calculations include both intramode (single mode) anharmonicity and intermode (mode coupling) rovibrational anharmonicity. The former is well understood, with, e.g., Morse-like “positive” anharmonicity acting to increase the number of states relative to the harmonic representation and quartic (or square-well) “negative” anharmonicity having the reverse effect. The effect of mode coupling is less well understood, although one may readily anticipate a variety of significant couplings, such as the coupling of torsions (internal rotations) with overall rotation and the coupling of stretches with the bends they make up. The quantitative effects of these couplings,



Page 4 of 50

4 including the cancelation of anharmonic effects with different signs, are not well understood, particularly at elevated temperatures. The present calculations allow for a detailed study of the quantitative effects of mode couplings and single mode anharmonicities, albeit classically. In practical kinetics and thermochemistry applications, MCPSI partition functions may be used semiclassically via Pitzer and Gwinn’s49 “useful approximation,” where the accurate (i.e., full-dimensional anharmonic rovibrational) partition function is approximated as the full-dimensional classical result corrected via the ratio of quantal and classical rigid rotor (RR) harmonic oscillator (HO) partition functions. Equivalently, the ratio of the exact classical result to the classical RRHO result can be interpreted as a semiclassical correction to the quantum mechanical RRHO model; this correction is the focus of the present study. Other semiclassical corrections, 69– 74 as well as quantal corrections via path integrals, could be employed in the context of MCPSI, but these are not considered here. This paper is organized as follows: In Sec. 2, our implementation of the MCPSI approach67,68 is reviewed, and modifications made here to facilitate its application to fluxional species are discussed. The analytic potential energy surfaces taken from Ref. 5 are also briefly reviewed. In Sec. 3, we report anharmonic correction factors, both with and without rotations, for a few representative species from five classes of systems: nonfluxional species, species with inversions but without torsions, species with methyl rotors, species with torsions and inversions, and species with secondary higher-energy minima. Qualitative features of the computed anharmonicity corrections are rationalized, and average anharmonicity corrections per mode at high temperatures are quantified. Section 4 is a summary.

2. Theory 2.1. MCPSI Our implementation of MCPSI has been described previously.67,68 Briefly, Monte Carlo coordinate sampling is performed using either Cartesian normal mode displacements or curvilinear “z-matrix” (stretch, bend, torsion) coordinates. The latter choice is particularly convenient when considering species with multiple minima, such as


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


5 several of the fluxional species considered here. When using curvilinear coordinates, the partial derivatives comprising the Jacobian are evaluated numerically, and in both cases the integration over momenta is performed analytically. Computed state counts using this approach were shown to be insensitive to the choice of coordinate system, with curvilinear coordinate sampling requiring up to 10x fewer sampled geometries than Cartesian normal mode displacement sampling for a given relative statistical uncertainty.67 The present set of results was obtained principally using curvilinear coordinate sampling. In fact, for many of the species considered here, sampling using Cartesian normal mode displacements about a single minimum would be prohibitively inefficient and/or require complex geometrical constraints to restrict partial sampling of equivalent adjacent minima. Insensitivity of the results to the choice of coordinate system was confirmed for systems where sampling using Cartesian normal mode displacements was feasible. In our previous MCPSI studies,67,68 a simple bias was used to sample geometries. In this approach, Nsamp geometries were selected from a hyperrectangle defined by singlecoordinate turning points and the energy ETP, where ETP is an input parameter. This hyperrectangle was divided into 100 nested hyperrectangles based on uniform divisions of each coordinate. The resulting nested shells were each sampled Nsamp/100 times, effectively biasing the sampled distribution as ~r instead of ~r!–3, where r is a generalized hyperradius of the sampled coordinates and !–3 is the number of vibrational degrees of freedom. This approach was modified slightly here to treat systems with multiple minima. For these systems, the nested hyperrectangle bias was removed for the one (or sometimes two) curvilinear coordinate(s) connecting the minima, e.g., for the torsional coordinate in CH3CH3. Biased sampling over the other coordinates was performed as before, with the removed coordinate(s) sampled without a bias over the allowed values. Without such a modification, a single minimum would be sampled more heavily than the others, degrading statistical convergence. We emphasize that using either strategy would eventually converge to the exact classical result with respect to Nsamp, albeit with varying efficiencies.



Page 6 of 50

6 As described previously,67 in our approach the rovibrational state count W is evaluated for a range of internal energies (0 – Emax) from the single ensemble of sampled geometries. Thermal partition functions are then computed by integrating dW/dE exp(!E/kBT) from 0 to Emax, where kB is Boltzmann’s constant. Often W can be converged to higher energies than the maximum sampled single-mode-displacement energy ETP, and Emax is therefore treated as an independent input parameter. Total sample sizes of Nsamp = 105 – 109 were used here, and this and other numerical parameters were chosen to converge the Q to better than 1%, unless otherwise indicated. In the largest production runs, Q from 200–3000 K was computed using 2160 cores running in parallel for 2 hours on Argonne’s Leadership Computing Resource Center’s Bebop cluster. For all the systems studied here, convergence tests with respect to the various numerical parameters were carried out, with an example of such a test is given in Sec. 3.1. More efficient and complex sampling strategies could certainly be designed to lower the required value of Nsamp. We do not pursue such strategies here, however, as even 109 analytic potential energy surface (PES) evaluations was not found to be prohibitively expensive. We define three thermal anharmonic correction factors:

f ! !1 Q / QHOQRR

(1)

fvib ! !1 Qvib / QHO

(2)

frot ! f / fvib ,

(3)

where " is the number of classically distinguishable lowest-energy minima connected via energetically-accessible paths, Q is the full-dimensional rovibrational partition function, Qvib is the rotationless anharmonic vibrational partition function, QHO is the harmonic oscillator partition function (for a single local minimum), and QRR is the asymmetric top rigid rotor partition function for the equilibrium geometry. Q, Qvib, QHO, and QRR are all evaluated classically, with Q and Qvib calculated using MCPSI and QHO and QRR calculated analytically. The rovibrational (total) correction f is of direct relevance to thermal kinetics and thermochemistry. The approximate separation of f into vibrational and rotational parts, as in Eqs. 1–3, is done to help rationalize the results


Page 7 of 50 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 presented below and to facilitate connections with other theoretical approaches where rotations are often neglected. As defined here, the correction factors do not include the significant errors that may arise due to anharmonicity in the zero point energy, and instead by construction f " 1 as the temperature T " 0. At finite temperature, one may anticipate that f scales, at least approximately, with the number of rovibrational degrees of freedom, ! = 3n – 3, where n is the number of atoms. We therefore consider the suitability of representing the present calculated anharmonicity factors as ! " 23

f = [1+ !(T )]

,

(4)

where # is the average relative anharmonic correction per mode at the temperature T, and we have assumed vibrations and overall rotations to contribute 1 and ! to the exponential scaling, respectively. No linear molecules are considered here, and so there are always three modes describing overall rotation. Corrections to Qvib can be used to isolate the vibrational anharmonicity, fvib = [1+ !vib (T )]! "3 ,

(5)

with the remaining rotational and rovibrational coupling anharmonicity expressed 3

frot = [1+ !rot (T )] 2 .

(6)

As demonstrated below in Sec. 3, the anharmonic corrections f, fvib, and frot show near linear dependence on T for some systems, while for other systems, #, #vib, and #rot are more complicated functions of temperature. Here, we use Eqs. 4–6 to quantify average per mode anharmonicities for all 22 species considered at consistent reference temperatures T to allow for systematic comparisons of the different molecular classes.

2.2. Potential energy surfaces The analytic PESs used here were published previously,5 where they were used to obtain high-accuracy anharmonic zero point energies from diffusion Monte Carlo calculations. Briefly, a Metropolis sampling algorithm was employed to generate a set of relevant ab initio energies for each system. These were least-squares fit using permutationally invariant polynomials 75 constructed from all atom–atom distances. Typically, the CCSD(T)/cc-pVTZ method was used, and we will also briefly consider fits



Page 8 of 50

8 obtained using larger and smaller basis sets. Energies of up to 50 kcal/mol (~17500 cm–1) above the classical minimum were included in the fits, with weights added to emphasize lower-energy regions. For a typical molecule, methanol, the surface was sampled with 50,000 ab initio energies. Near the minimum, these data were fit with errors of just 1 cm–1, with an overall fitting error of 20 cm–1. A strategy was implemented for constraining sampling of the fitted PESs to physically relevant regions. Unlike in the previous DMC study, here we often sample the fitted PESs at geometries with large coordinate displacements in more than one degree of freedom. These geometries can be deep in the repulsive wall and outside the region for which the PESs were originally parameterized. At these geometries, the fitted PESs sometimes returned erroneously low energies instead of correctly returning high energies. Therefore, for each sampled geometry, a test energy Vtest was evaluated. If the test energy was too high, the value of the analytic PES was discarded, and instead the contribution to W was set to zero (simulating the result of a high-energy V). The criterion for rejection used here was Vtest > max(EX, mV), where EX and m are input parameters. A convergence test with respect to these parameters is presented in Sec. 3.1. For systems with a single minimum, Vtest was chosen to be the usual harmonic expansion of the fitted potential energy surface about the minimum. For systems with multiple minima, the coordinate connecting the minima (e.g., the torsion) was left out when evaluating the otherwise harmonic Vtest.

3. Results and Discussion 3.1. Nonfluxional systems: H2O, HO2, CH2O, CH4, and CH2CH2 We first consider a set of five systems without nonlocal torsions and inversions and with sizes ranging from 3 to 6 atoms and 6 to 15 rovibrational degrees of freedom. The vibrational and rovibrational MCPSI anharmonic correction factors are shown in Fig. 1. Results for the two corrections are qualitatively similar to one another, both scaling linearly with T and exponentially with the number of degrees of freedom, as in Eqs. 4–6. Despite the different atomic compositions and different relative numbers of bends and stretches, these five systems show similar average anharmonicities per mode. The curves in Fig. 1 are well described by Eqs. 4 and 5 with # = 2.7–3.5% and #vib = 2.7–3.3% for T


Page 9 of 50 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


9 = 2500 K. Using Eq. 6, we can similarly quantify the contribution to the total anharmonic correction f from the rotational modes. For the five systems in Fig. 1, the rotational anharmonicity per mode is approximately ! to " as large as for the vibrational modes, on average. This result was anticipated in our Eqs. 4–6 where the rotational degrees of freedom were each assigned exponents of !. The resulting values of #rot = 2.8–4.3% (T = 2500 K) therefore have similar magnitudes as #vib. The per mode anharmonicities for these systems, and all the systems considered in this work, are summarized in Table 1 for two temperatures. Before moving on to fluxional systems, we demonstrate numerical convergence and comment on the expected overall accuracy of the present approach. In Fig. 2, fvib for CH2O is shown for various choices of the numerical parameters ETP and EX using curvilinear coordinate sampling. For ETP = 15000 cm–1 and when no harmonic test potential is used to constrain sampling of the fitted PES to physically relevant regions (i.e., when EX = #), the calculated value of fvib varies from 13–1.5 from 200–2500 K and cannot be seen in the figure. This unphysical behavior results from sampling highly attractive regions of the fitted PES outside the region where it was parameterized. This error is eliminated by choosing a smaller sampling domain (i.e., by lowering ETP) or by implementing the harmonic test potential described in Sec. 2.2 (here, with m = 0 and EX = ETP or 2 ETP). For these choices, the computed values of fvib all agree fairly well with one another at low T, and for some choices (e.g., ETP = 15000 cm–1, m = 0, and EX = 2 ETP) convergence is demonstrated up to 2500 K. Similar tests were done for every system considered in this article, and often this same set of numerical parameters was found to be appropriate. For some systems, an energy-dependent test potential screening criterion was found to be necessary, with ETP = 15000 cm–1, m = 2, and EX = 10000 cm–1. The CH2O system is one of the few systems considered here where accurate results can be conveniently obtained by sampling in Cartesian normal mode displacements as well as curvilinear coordinates. Figure 2 includes a demonstration of the independence of the present results with respect to this choice, with the curvilinear and Cartesian normal mode displacement results agreeing with one another to better than 0.5%.



Page 10 of 50

10 Finally, we note that all the curves shown in Fig. 2 were calculated using a relatively small number of sampled geometries (Nsamp = 105). The small deviations of the results for the several “converged” calculations in Fig. 2 are on the order of 0.5%, which may be assigned as the statistical uncertainty for this system and this choice of Nsamp. Elsewhere in this work Nsamp = 106!109. The results for each system were confirmed to differ from those for Nsamp/10 by less than 2% suggesting statistical error bars of less than 1% for the results reported in Fig. 1 and below. Next, we briefly consider the absolute accuracy of the present approach. In our previous study,67 MCPSI results for fvib for CH4 were compared with vibrational configuration interaction (VCI) results calculated using the same PES. The two results differed at low T, with the MCPSI result showing a larger anharmonic correction than the quantum mechanical result. At 300 K, for example, fvib = 1.031 and 1.002 for the MCPSI and VCI methods, respectively. While this is a large relative error (a factor of 15), it is a small absolute one (just 3% in the partition function). Both the relative and absolute errors in the MCPSI approach were shown to decrease with increasing temperature, with the two approaches agreeing quantitatively at 1500 K (fvib = 1.191 and 1.200 for the MCPSI and VCI methods, respectively). In another test of the Pitzer–Gwinn approach for a full dimensional polyatomic system at high temperatures, Lynch et al.76 compared rovibrational partition functions for H2O2 computed using their (exact) path integral approach with its classical limit. When the effect of anharmonicity in the zero point energy is removed, their classical limit is equivalent to the present approach. They found differences in f (again excluding zero point effects) of 7% at 300 K and just 2% at 2400 K. These results are likely general, at least qualitatively, where the classical approach may be expected to overpredict quantum mechanical results at low temperatures, when f is typically close to unity, and to be more accurate at higher temperatures where f is typically largest. Note also that for species with H atoms connected to a central heavy atom, the quantum correction to the classical rotational partition function at low temperatures may have a magnitude similar to the rovibrational corrections being quantified here.77 Finally, we consider basis set dependence of the calculated partition functions. For the nonfluxional species in this subsection, the harmonic frequencies show only a


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


11 weak dependence on basis set. For CH2O, for example, QHO for the cc-pVDZ and ccpVTZ basis sets differs from the cc-pVQZ result by 5% and 2.4% respectively. Similar relative differences were found for classical anharmonic partition functions Q computed using fitted analytic PESs parameterized against CCSD(T) energies calculated using these three basis sets. The ratio f therefore benefits significantly from error cancellation, with errors in f for the smaller basis sets reduced by nearly a factor of 10. Specifically, the calculated values of f for the cc-pVDZ and cc-pVTZ basis sets differ from the ccpVQZ results by less than 0.6% and less than 0.3%, respectively, over the temperature range considered. This suggests the practical strategy of computing anharmonic corrections using a lower level of theory than may be required to converge the partition functions individually. A similar conclusion regarding strategies for anharmonic zero point energy calculations was given in Ref. 5, where particularly slow convergence in the harmonic frequencies with respect to basis set was observed for some fluxional species. 3.2. Fluxional systems with inversions: 3CH2, CH3, NH3, and CH2CH Anharmonic corrections for four systems with inversions but without torsions are shown in Fig. 3. For this set of systems, CH3 is a clear outlier due to the negative anharmonicity associated with the strong quartic component of its umbrella motion.78,79 As may be anticipated, the negative anharmonicity for this lowest frequency mode results in f < 1 at low temperature. At higher temperatures, the positive anharmonicities associated with the remaining five bends and stretches somewhat cancel out the effect of the umbrella mode, such that fvib and f deviate from unity by less than 15% and 10%, respectively, over the range of temperatures shown. For the other three systems, the anharmonicity corrections again increase with the number of degrees of freedom (as noted for the nonfluxional species in Sec. 3.1) but with a T-dependence that is not as linear as observed in Fig. 1. The per-mode vibrational anharmonicity corrections are also larger than those in Fig. 1, with #vib = 3.8–5.3% at 2500 K. Both the larger per-mode anharmonicities and the observed nonlinearity with T may be associated with the presence of mirror image wells (" = 2) with inversion barriers of 1826, 2069, and 2238 cm–1, and for CH2CH, CH2, and NH3, respectively. Critical point



Page 12 of 50

12 energies for all of the species considered here are collected in Table 2. These barriers are accessible at elevated temperatures, as demonstrated next for CH2. Similar results were found for CH2CH and NH3 but are not shown. Note that HO2 in the previous section also has mirror image wells, but with a much higher inversion barrier (~40 kcal/mol). An additional calculation was carried out where only the curvilinear bending coordinate of CH2 was sampled, with the remaining degrees of freedom frozen at their equilibrium values. The computed anharmonicity corrections to the full-dimensional vibrational and one-dimensional bending densities of states #vib are shown in Fig. 4. With the bending coordinate isolated, the effect of the barrier can be clearly seen as the density of states rises sharply at energies just below the classical threshold energy and falls sharply at energies just above it. This behavior near a classical threshold has been analyzed in detail in the context of hindered rotors.52,54 Below threshold, the one-dimensional bending anharmonicity corrections exhibit Morse-like positive anharmonicity due to the flattening out of the potential for near-linear bond angles (or near-planar umbrella angles in the case of NH3). Above threshold, the systems exhibit negative anharmonicity, which may be attributed to two sources. First, the double well potentials feature a significant quartic component. Second, the analytic harmonic reference density of states #HO assumes integration over the bending coordinate from –$ to $, whereas the range of allowed values is actually 0 to 180o. The effect of removing the contribution from disallowed values of the bending angle from #HO is shown in Fig. 4. In one dimension, the onset of this effect occurs around 6500 cm–1, which is the energy of the harmonic bending potential for an angle of 180o. Notably, with this contribution removed, the resulting state density appears nearly harmonic (independent of energy) above 6500 cm–1, with an effective harmonic frequency for half of the double well that is ~20% lower than the single well frequency. Results for the full-dimensional vibrational density of states are also shown in Fig. 4, where the two stretches are included along with the bend. The effect of the bending coordinate on the full-dimensional #vib gives rise to similar qualitative features in the three-dimensional system, although the observed features are shifted to larger total energies, are less sharp, and can no longer be quantitatively associated with the inversion barrier height. Convolving the one-dimensional anharmonic bending density of states


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


13 with harmonic densities of states for the two stretching modes reproduces this shape qualitatively but with an anharmonic correction of approximately " to $ the magnitude. The mitigation of sharp single-mode classical features in higher-dimensional density of states convolutions has been discussion previously, again in the context of hindered rotors.80,81 Next we consider rovibrational anharmonic corrections to the partition function, as shown in Fig. 3(b). Again, CH3 is an outlier. Notably, when rotations are included for the other species, the computed anharmonicity corrections become more similar to those for the nonfluxional species quantified in Sec. 3.1. Specifically, the rovibrational anharmonic corrections for CH2, CH2CH, and NH3 vary more linearly with T than the vibrational corrections, with per-mode anharmonicity corrections of # = 2.6–4.0% for T = 2500 K, which may be compared with # = 2.7–3.5% for the nonfluxional systems in Fig. 1(b). This result suggests that the qualitative differences in fvib noted for nonfluxional systems in Figs. 1(a) and the systems with inversions in Fig. 3(a) are due, at least in part, to the separation of vibrations and rotations. Quantitatively, as shown in Table 1, fluxional systems with inversions have only slightly larger overall anharmonicities per mode than nonfluxional ones, despite the much larger vibrational contribution per mode and much lower rotational contribution. This may be rationalized by noting that systems with inversions typically feature qualitative geometry changes as their inversion barriers are surmounted, which can introduce significant rovibrational coupling. An extreme case of course is a system with access to linear geometries (such as CH2), where a rotation is converted into a vibration as collinearity is approached. When rovibrational couplings are fully included, the per-mode anharmonicity corrections for species with inversions appear similar to those for nonfluxional species. The reduction in the magnitude of the anharmonic effect of the bending mode can be seen in the rovibrational state densities shown in the inset in Fig. 4 for CH2. With overall rotations included, the rovibrational density of states for the bend appears qualitatively different than #vib, and in fact it appears approximately harmonic (independent of energy) at low energies. This is an example of the general difficulty in constructing accurate full-dimensional rovibrational properties from separable one-



Page 14 of 50

14 dimensional ones as mode couplings can introduce qualitative changes, such as the strong coupling of a bend with overall rotation that has been observed here. The use of curvilinear coordinates greatly aids in the present applications. As noted previously,67 convergence of the anharmonic corrections for the vinyl radical required ~10x fewer sampled geometries when using curvilinear coordinates than when sampling using Cartesian normal mode displacements. Furthermore, when using Cartesian sampling, one of the two equivalent minima must be chosen as the reference geometry, and a geometric constraint must be included to avoid sampling the other well without any loss of symmetry. The implementation of such a geometric constraint is not always obvious for many nonfluxional systems. For example, we tested a system-specific cutoff based on an effective umbrella angle for NH3, but the results showed seemingly systematic deviations from the curvilinear result. This comparison was made difficult by the poor sampling efficiency using this approach, and a more detailed study was not pursued. System-specific geometric cutoffs are not needed when sampling using curvilinear sampling. As shown in Table 1, the “per mode” vibrational anharmonicity corrections for the systems quantified so far vary from #vib = 1.4–3.9 at 1500 K and 2.7–5.3 at 2500 K, excluding CH3, where #vib % –2 at elevated temperatures. Physically, these corrections result from the typically positive (Morse-like) single mode anharmonicity associated with each stretch and bend, as well as the positive or negative anharmonic couplings between modes. A simple approach sometimes employed to approximate anharmonic effects is to use the quantum mechanical harmonic oscillator formulas but with fundamental frequencies used in place of the usual harmonic ones.76 By construction, this approach accurately builds in anharmonicities at the low energies relevant to the 0 " 1 transitions so long as the fundamentals are accurately known. We quantified how much of the anharmonicity corrections this simple approach recovers at 1500 and 2500 K, with fundamental frequencies computed using the VPT2 calculations reported in Ref. 5. This approach corresponds to an additional multiplicative factor in the reference partition QM QM / QHO function (i.e., in the denominators of Eqs. 1 and 2) of QHO(fund) , where the harmonic

oscillator partition functions are evaluated quantum mechanically relative to the zero point energy using fundamental and harmonic frequencies, respectively. By construction


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


15 this ratio equals unity at 0 K, and, because fundamental frequencies are typically smaller than harmonic frequencies due to the predominance of positive Morse-like anharmonicity, this ratio typically increases with T. CH3 is again an exception due to the strong negative anharmonicity associated with its quartic umbrella potential. Quantitatively, this approach recovers roughly half (43–63%) of the total vibrational anharmonicity correction, with similar results at 1500 and 2500 K, again excluding CH3, for which the approach fails qualitatively predicting positive instead of negative anharmonicity corrections. Lynch et al. came to similar conclusions comparing this approach with path integral results.76

3.3. Systems with methyl rotors: CH3CH3, CH3OO, CH3OH, CH3CHO, CH3CO, and 3

CH3CH Calculated values of f and fvib for six systems with symmetric three-fold periodic

methyl rotors (" = 3) are shown in Fig. 5. These systems have 6 to 8 atoms with 15 to 21 rovibrational degrees of freedom, and they fall into three groups based on their torsional barriers. Ethane is the lone example of a strongly hindered system, with a barrier height (980 cm–1) large enough to significantly hinder free rotation even at elevated temperatures. The inflection point around 700 K indicates the temperature at which free rotation begins to influence f and fvib. This free rotation contribution has some effect on lowering f and fvib above 700 K, but f continues to increase with T even up to 2500 K. In fact, the per-mode anharmonicities for ethane are not much smaller than those quantified for the systems in Secs. 3.1 and 3.2, with # = 2.6% at 2500 K, suggesting a relatively small influence of the free-rotation torsional regime for this system. In the other limit, the systems CH3CO and CH3CH are only very weakly hindered with low torsional barriers of 161 and 115 cm–1, respectively. The sharp turnover in f and fvib at temperatures below 100 K can be associated with ready access to free rotation, even at the low energies populated at 100 K. Free rotation in one-dimension has an analytic T–1/2 temperature dependence, and similar behavior can be clearly seen in the present full-dimensional results. Notably, the negative anharmonicity from this single mode dominates the positive anharmonicities associated with the remaining rovibrational



Page 16 of 50

16 modes, even at high temperatures, giving rise to a negative temperature dependence overall. The third group consists of CH3OO, CH3OH, and CH3CHO, which all have similar torsional barriers of 336, 384, and 403 cm–1, respectively. Near room temperature and below, the majority of highly populated energies are below the threshold for internal rotation, and the anharmonicity correction increases sharply with T as more and more of the flattening potential near threshold is accessible. At higher energies, increased access to free rotation introduces a negative temperature dependence, which is seen here to nearly quantitatively cancel out the positive anharmonicities in this mode at lower energies and associated with the remaining modes. The result of this cancellation is that the fully harmonic partition function for these systems deviates from one-third (" = 3) of the fully anharmonic vibrational partition function by less than ~15% for temperatures above 1000 K. A simple one-dimensional periodic rotor is used next to rationalize the trends in Fig. 5. We assume a hindering potential of VHR($) = Vb/2 (1 – cos "$) for the torsion angle $ and barrier higher Vb, and the one-dimensional hindered rotor partition function is computed numerically:52,54 !

QHR = kBT / ! A " d" exp(!VHR / kBT ) . 0

(7)

As is frequently done, the effective moment of inertia A is defined as the reduced moment of inertia of the two rotating groups at their equilibrium geometries about the axis of their connecting bond.82 Differentiating the hindering potential allows one to calculate an effective harmonic frequency for the torsion (equal to ! AVb ), from which an effective eff function may be calculated. By construction, the ratio harmonic partition QHO eff fHR = QHR / QHO

(8)

tends to unity at low temperature. This ratio is a hindered rotor correction to the harmonic reference partition functions (i.e., the denominators) in Eqs. 1 and 2. Equation 8 was evaluated using the barrier heights given in Table 2 and CCSD(T)/cc-pVTZ geometries for A, and the results are shown in Fig. 6(a). The onedimensional corrections fHR are seen to qualitatively reproduce several trends observed in the full-dimensional curves in Fig. 5, including the temperatures at which the transition


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


17 from positive to negative anharmonicity occurs, as well as the relative magnitudes of the two effects. In Fig. 6(b), the ratio f/fHR = Q/QHOQRRfHR is shown. The average correction per mode to the reference partition function QHOQRRfHR can be quantified using an equation similar to Eq. 4, giving # = 2.3–3.4% at 2500 K. We therefore conclude that once the effect of the hindered rotor is accounted for via Eq. 8, the remaining anharmonicities for six systems are fairly similar to one another and to those quantified for the systems without torsions in Secs 3.1 and 3.2. The ratio plotted in Fig. 6(b) is relevant to many practical applications, where torsions are sometimes treated as one-dimensional hindered rotors while other modes are treated harmonically. For all six systems, employing a hindering potential improves the accuracy of the reference partition function at low temperatures. At higher temperatures, treating the torsion for ethane as a hindered rotor improves the reference partition function (lowering the magnitude of the correction in Fig. 6(b) relative to its magnitude in Fig. 5(b)), although the improvement is small and decreases with temperature. For the three systems shown here with intermediate torsional barrier heights, however, the fully harmonic model is found to be more accurate than the hindered rotor model at temperatures above ~1000 K and is 3–4x more accurate at 2500 K. For the two systems with nearly free rotors, the two models have errors with similar magnitudes but opposite directions relative to unity. The counterintuitive result that systems like CH3OH are more accurately described at high temperatures without making a hindered rotor correction relies on a cancellation of errors––specifically a cancellation of the negative high-temperature torsional anharmonicity and the positive stretching and bending anharmonicities for the other modes. As seen here for the other methyl rotor systems in this section, however, this cancellation cannot be relied upon, and it may be difficult to anticipate the overall effect quantitatively. Even in this small sample of species with common organic chemical structures, torsional barriers vary enough so as to span the nearly harmonic and nearly free rotor regimes. For species in the intermediate regime, the relative importance of free and hindered rotation is temperature dependent around 300 K. This result, along with the significant rovibrational coupling at high temperatures quantified in Sec. 3.2, highlights



Page 18 of 50

18 the general difficulty in extrapolating room temperature results (at which the majority of experimental determinations are available) to higher temperatures. Furthermore, these results demonstrate that two seemingly similar systems, such as CH3CO and CH3OO, may show qualitatively different behavior in their rovibrational properties. We note that the average per mode rotational anharmonicity corrections for these torsional systems (#rot = 3.4 – 5.2%) are just slightly larger than those quantified for the nonfluxional systems in Sec. 3.1. While one may generally anticipate larger rovibrational coupling for systems with internal rotors, the presence of symmetric methyl rotors evidently mitigates the effect here. CH3CH is an example of a class of fluxional systems where the torsional mode is strongly coupled to an energetically accessible inversion of one of the rotating groups. Both motions can result in conformer changes, suggesting a two-dimensional representation for describing this periodic system instead of the usual one-dimensional hindered rotor representation. More complex examples of this situation are considered in the next section. Here, we present a detailed study of the coupled inversion and torsional modes for CH3CH, which has a fairly high barrier to inversion (2660 cm–1) thus simplifying the analysis. A two-dimensional contour plot of the fitted anharmonic PES as a function of the C–C–H bending angle and the H–C–C–H torsional angle is given in Fig. 7(a), along with a harmonic expansion about one of the three equivalent minima. The PES in both coordinates features significant anharmonicity. The torsion is clearly more accurately described as a periodic hindered rotor, while the bend appears “Morse-like” with a steeper-than-harmonic hard wall and a flattening potential as the bend approaches linearity. As plotted in Fig. 7(a), the two anharmonic coordinates appear to be largely separable from one another. Because the harmonic model is a poor description of the potential dependence of both modes, one may be tempted to propose an improved separable model as the sum of a hindered rotor and Morse-like potential. Implicit in any separable model is the assumed separability of the kinetic energy as well as the potential energy. The torsion and inversion coordinates considered here are strongly coupled kinetically. For smaller and smaller deviations from linearity, a given change in the torsional angle corresponds to smaller and smaller geometry changes, and


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


19 at collinearity the torsion cannot carry any kinetic energy. This is equivalent to the observation that the slice through the contour plot in Fig. 7(a) for collinear C–C–H corresponds to just a single geometry. The neglect of this coupling (and of the kinetic coupling for any two curvilinear modes) can be a significant source of error, just as significant as the neglect of potential energy coupling. The non-rectilinear nature of these coordinates and their kinetic coupling is represented by the Jacobian, which describes the mapping of curvilinear coordinates onto (kinetically separable) rectilinear ones. The Jacobian Ji appears as a weight in the Monte Carlo sampled integrand, e.g.,

Q = C1 " J(xi )exp(!V (xi ) / kBT ) , i

(9)

where i labels the sampled geometry xi, and C is a normalization constant. One may define an effective temperature-dependent potential for which the system is kinetically separable by bringing the Jacobian into the exponential in Eq. 9, and writing

Q = C1 # exp(!V "(xi ) / kBT ) ,

(10)

V ! = V " kBT ln J ,

(11)

i

where

and the corresponding effective Jacobian, J', is unity. The relative magnitude of J as a function of the C–C–H angle is plotted in Fig. 7(b), where it is shown (necessarily) to decrease sharply to zero for angles approaching linearity. The large variation in J reflects the large kinetic coupling for these two modes. Figure 7(b) also compares the anharmonic bending potential, its harmonic expansion, and the kinetically separable effective potentials defined by Eq. 11 for three temperatures. All of the models were zeroed to agree at the equilibrium geometry. The kinetic separability correction qualitatively changes the shape of the potential. Notably, the resulting effective potential appears to agree more closely with the harmonic potential at high temperatures that with the anharmonic potential. This suggests that treating the torsion as a hindered rotor and the bend as a harmonic oscillator may be an accurate separable curvilinear representation of these two modes, at least at high temperatures. “Improving” the harmonic representation of the bend to look Morse-like could in fact degrade the accuracy of the model if kinetic coupling is also neglected. This is likely a general result



Page 20 of 50

20 suggesting caution when using curvilinear coordinates in reduced dimensional models. While curvilinear coordinates often greatly reduce potential couplings, the corresponding increases in kinetic couplings, which are more difficult to quantify, can be significant.

3.4. Fluxional systems with asymmetric torsional barriers and inversions: H2O2, NH2NH2, CH2OH, and CH3CH2 In this section, we consider systems with more complex vibrational structures than discussed in previous sections. We therefore do not attempt detailed quantitative rationalizations of the results for these systems. Instead we quantify the average per mode anharmonicities and include brief qualitative comments. The H2O2 and NH2NH2 systems both feature torsions with asymmetric barriers to rotation. The higher-energy barriers are similar in the two systems (~2750 cm–1) and likely prevent significant contributions to the partition function from fully free rotations. The lower-energy barriers are accessible, however, with the barrier for H2O2 (365 cm–1) lower than that for NH2NH2 (748 cm–1). As shown in Fig. 8, the magnitude and temperature dependence of the anharmonic corrections for H2O2 are quantitatively very similar to the methyl rotor systems with barriers around 350 cm–1. This similarity is despite the different rotational symmetry " and high-energy barrier to free rotation in H2O2. The fully harmonic model is very accurate for H2O2, again presumably due to the cancellation of errors described in Sec. 3.3. The anharmonicity corrections for NH2NH2 are the largest quantified here so far. Along with the low-frequency torsion, this system may also undergo –NH2 inversions analogous to the umbrella motion in NH3 and with a similar barrier height (2247 cm–1). A plot of the torsional mode and one of the umbrella modes is given in Fig. 9. This system of three coupled fluxional modes gives " = 8 and could explain the relatively large anharmonicity for this system. We also note that NH3 showed larger anharmonicities than the other non-torsional species in Secs. 3.1 and 3.2, which might indicate that N-centered stretches and bends generally have larger anharmonicities than C-centered ones. Next, we consider two species that feature low-energy inversions, as highlighted previously in the anharmonic zero point energy study.5 As discussed there, CH2OH is nearly planar with a low-energy barrier (139 cm–1) for inversion and a larger torsional


Page 21 of 50 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


21 barrier (1670 cm–1) separating " = 4 classically distinguishable minima. For CH3CH2, " = 6 classically distinguishable minima are separated by a mixed inversion/torsional barrier with a very low energy of just 32 cm–1. These values of " may not be evident from the chemical structures alone, and potential energy slices demonstrating these symmetries were presented in Ref. 5. Notably, for both systems, the inversion barriers are much lower than in the CH3CH system analyzed in detail above, further complicating the vibrational structure of the coupled inversion and torsion. The anharmonicity correction for CH3CH2 shows a strong negative temperature dependence, as shown in Fig. 8, lowering the partition function relative to the harmonic model for the entire temperature range considered here and by as much as " at 2500 K. This effect is approximately twice as large as observed for the two nearly free methyl rotor systems considered in Sec. 3.3. This suggests that, except at very low temperatures, the CH3CH2 system behaves similar to the CH3CH and CH3CO systems, with an effective rotational symmetry number of 3. For CH2OH, the anharmonic correction is nearly negligible over the entire temperature range considered here. The potential energy describing the inversion motion has a significant quartic component, contributing a negative anharmonicity for this lowest-frequency mode that largely cancels the positive anharmonicities associated with the remaining bends, stretches, and torsion. This is similar to the situation for CH3, with the larger number of modes for CH2OH likely responsible for the more complete cancellation found here than for CH3. The rotational anharmonicity corrections per mode for these four systems vary from #rot = 4.3–5.9% at T = 2500 K, which is similar to #rot for the methyl rotor systems in Sec. 3.3 and somewhat larger than nonfluxional systems in Sec. 3.1. As noted in Sec. 3.2, systems with accessible inversions are likely to feature relatively larger rotational anharmonicities and couplings.

3.5. Systems with secondary minima: CH2CHOH, HC(O)OH, and NH2OH Finally, we consider three systems with asymmetric torsions and non-equivalent secondary minima. The most complex of these is NH2OH, having two fluxional modes, one rotation and one inversion. A plot of these two modes is given in Figure 10. The



Page 22 of 50

22 vibrational and rovibrational anharmonic correction factors for these systems (Fig. 11) are large, with average per mode anharmonic corrections of # = 8.7–9.0% for T = 2500 K. These large values arise principally from vibrational anharmonicity (#vib = 7.9–8.6%), whereas the rotational anharmonicities (#rot = 2.6–5.9%) are similar in magnitude to the systems with torsions considered above. These results may be largely explained by the presence of the higher-energy isomers. The secondary minima for CH2CHOH, NH2OH, and HC(O)OH have energies of 539, 1512, and 1565 cm–1, respectively. These are low enough to be populated at high temperatures, whereas the reference harmonic partition functions in Fig. 11 represent just the single (" = 1) low-energy minimum. The barriers to free rotation are fairly large for these species, and, although the conformers may be chemically distinct at low enough temperatures, we make no effort to compute conformer-dependent partition functions. Instead, we improve the reference harmonic partition function to include the secondary minimum, (2) Q2HO = QHO + Q (2) / kBT ) , HO exp(!V

(12)

(2) where Q (2) HO is QHO computed using the frequencies at the secondary minimum, and V

is its energy relative to the first minimum. We also consider the RRHO analog of Eq. 10 and the simpler model Q2HO ! QHO (1+ exp("V (2) / kBT )) ,

(13)

where frequency changes at the higher-energy minimum are neglected. Vibrational anharmonicity corrections for the reference partition functions in Eqs. 12 and 13 are shown in Fig. 12. Including the effect of the secondary minimum reduces the average per mode anharmonicity corrections to # = 3.9–5.1% for T = 2500 K. These values are still relatively high compared to the systems considered above, suggesting that representing the torsional mode as a sum of harmonic contributions is a source of significant error. One could employ an asymmetric hindered rotor correction instead of Eq. 12, but this was not attempted here. Comparing the results using Eqs. 12 and 13 (shown in Fig. 12) can be used to quantify one effect of vibrational coupling of the nontorsional modes to the torsional one. For CH2CHOH and NH2OH, the effect is small, while a fairly large effect is seen for HC(O)OH. This system features an internal


Page 23 of 50 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


23 hydrogen bond in its lower-energy geometry but not its higher-energy one, leading to significant frequency changes in the non-torsional modes for the two minima. A variety of methods have been developed for incorporating adiabatic frequency changes and asymmetric torsional potentials,45,58,74 and a comparison of these methods with the present full dimensional results could be made in a future study.

4. Conclusions One goal of this study was to demonstrate the usefulness of Monte Carlo phase space integrals (MCPSI) for characterizing fully anharmonic rovibrational properties at elevated temperatures and energies. The principal advantage of MCPSI is that it requires no a priori assumptions about the underlying vibrational structure, and it can therefore be applied to systems with difficult to characterize and strongly coupled vibrational modes. Twenty-two molecules and radicals were considered here, with up to 21 rovibrational degrees of freedom and with several systems featuring torsions and other low-frequency strongly-coupled fluxional modes. Trends in the results were quantified and qualitatively interpreted, with an emphasis on characterizing phenomena likely to be more significant at high temperatures than in the relatively more widely studied low temperature regime. For the nonfluxional systems included in the test set, the computed anharmonicity corrections at high temperatures were shown to be fairly consistent and independent of chemical composition, increasing the partition function by ~3% per mode at 2500 K, on average. A simple model where fundamentals are used in place of harmonic frequencies recovers half of this anharmonicity. For systems with inversions and without torsions, similar values were obtained, but only for the fully coupled rovibrational partition function. The separation of rotations and vibrations for these systems was shown to complicate the interpretation of vibrational anharmonicity at high temperature due to the relatively large rovibrational coupling associated with inversion-induced geometry changes. Several systems with methyl rotors were also studied, and their anharmonicity corrections were found to be a strong function of the torsional barrier height. These trends were normalized via a simple torsional correction, and the remaining modes showed a positive anharmonicity with magnitudes similar to those quantified for systems without torsions (again ~3% per mode at 2500 K).



Page 24 of 50

24 The consistency of these corrections confirms the empirical observation that replacing harmonic frequencies with fundamental frequencies (which typically exhibit positive anharmonicities) often improves the accuracy of harmonic oscillator predictions of rovibrational properties. However, for several systems studied here which feature modes with significant negative anharmonicity (such as quartic bends and torsions), the harmonic oscillator representation using harmonic frequencies was already found to be accurate due to the near quantitative cancellation of positive and negative anharmonicities. Models that incorporate only one kind of anharmonicity (such as an improved treatment of the torsional mode) can in fact be less accurate than fully harmonic models, as these cases demonstrate. Finally, several systems with strongly coupled low-frequency large amplitude modes were considered. These systems showed per-mode anharmonicity corrections that varied significantly from –7% to 9% per mode at 2500 K. Simple models for predicting these effects are not immediately evident. For example, several of these highly fluxional systems feature a torsion coupled strongly to one or more inversion modes of the rotating fragments. An analysis of separable anharmonic models for these modes was given to demonstrate the significant temperature dependence in kinetic couplings, suggesting that the accurate low temperature separable models may not be suitable at high temperatures. Disadvantages of the MCPSI approach include the neglect of quantum effects and challenges in sufficiently sampling the high-dimensional vibrational spaces. Progress regarding the former challenge was not made here, although the good accuracy of Pitzer and Gwinn’s semiclassical approximation adopted here has been demonstrated elsewhere at moderate and high temperatures.52,56,67,76, 83 Quantifying the accuracy of other semiclassical models and incorporating path integrals may be useful as future work. Regarding the latter challenge, we have demonstrated the convenience of sampling in curvilinear coordinates, enabled very generally in our approach via numerical evaluations of the required Jacobians. This strategy allows for the convenient treatment of systems with one or more low-frequency motions. A variety of sampling efficiency improvement strategies, although not needed in the present applications with as many as 21 degrees of freedom, could be implemented to treat larger systems. Alternatively, our own previously proposed strategy68 of representing the full-dimensional density of states via convolutions


Page 25 of 50 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


25 of reduced dimensional single mode, pairwise, 3-mode, etc., intrinsic state densities, which greatly reduces the computational cost, could be used.

Acknowledgments

This work was supported by the U. S. Department of Energy, Office of Basic Energy Sciences, Division of Chemical Sciences, Geosciences, and Biosciences through Argonne National Laboratory. Argonne is a U. S. Department of Energy laboratory managed by UChicago Argonne, LLC, under Contract Number DEAC02-06CH11357. Z. G. was supported by the U. S. Department of Energy Office Science’s Summer Undergraduate Laboratory Internship (SULI) program. Software development was supported by the Exascale Computing Project (ECP), Project Number: 17-SC-20-SC, a collaborative effort of two DOE organizations, the Office of Science and the National Nuclear Security Administration, responsible for the planning and preparation of a capable exascale ecosystem including software, applications, hardware, advanced system engineering, and early test bed platforms to support the Nation's exascale computing imperative. We gratefully acknowledge computing resources provided by Blues and Bebop, two high-performance computing clusters operated by the Laboratory Computing Resource Center at Argonne National Laboratory.



Page 26 of 50

26 Table 1. Average anharmonicity corrections per mode (%) at two temperaturesa System # #vib #rot Sec. 3.1 H2 O 1.4 [2.7] 1.4 [2.7] 1.5 [2.8] HO2 1.7 [3.3] 1.8 [3.3] 1.7 [3.3] CH2O 1.5 [3.0] 1.5 [2.9] 1.4 [3.2] CH4 1.7 [3.5] 1.7 [3.3] 2.0 [4.2] CH2CH2 1.5 [3.0] 1.4 [2.9] 2.0 [4.3] Average 1.6 [3.1] 1.6 [3.0] 1.7 [3.6] Sec. 3.2 CH3 –1.1 [–0.7] –2.0 [–2.2] 2.5 [5.7] 3 CH2 1.6 [2.6] 3.9 [5.3] –2.8 [–2.5] NH3 2.6 [4.0] 2.8 [4.3] 1.6 [3.0] CH2CH 2.2 [3.8] 2.3 [3.8] 1.4 [3.4] b Average 2.1 [3.5] 3.0 [4.5] 0.0 [1.3] Sec. 3.3 CH3CH3 1.9 [2.6]; 1.1 [2.5]c 1.8 [2.4] 2.4 [5.2] c CH3CHO 0.8 [0.8]; 0.8 [2.3] 0.7 [0.5] 1.8 [3.7] c CH3OH 1.2 [0.8]; 1.3 [2.9] 1.0 [0.3] 2.2 [4.8] c CH3OO 1.1 [1.0]; 1.6 [3.4] 1.0 [0.5] 2.3 [4.9] c CH3CO –1.7 [–2.2]; 0.7 [2.6] –2.2 [–3.0] 2.2 [4.8] 3 c CH3CH –2.4 [–3.6]; 1.1 [2.3] –2.8 [–4.4] 1.6 [3.4] c Average 0.1 [–0.1]; 1.2 [2.9] –0.1 [–0.6] 2.1 [4.5] Sec. 3.4 H2 O2 1.9 [1.5] 1.8 [0.7] 2.5 [4.7] NH2NH2 4.5 [6.1] 4.6 [6.1] 3.1 [5.9] CH2OH 0.2 [–0.5] –0.1 [–1.3] 2.0 [4.3] CH3CH2 –6.0 [–7.2] –6.8 [–7.3] 2.1 [4.4] Sec. 3.5 CH2CHOH 6.0 [8.7]; 2.9 [4.7]d 5.5 [8.2] 0.1 [4.2] d HC(O)OH 3.8 [8.8]; 1.8 [3.9] 3.4 [7.9] 1.3 [2.6] d NH2OH 4.7 [9.0]; 3.1 [5.1] 4.4 [8.6] 2.8 [5.9] d Average 4.4 [8.2]; 2.6 [4.5] 4.8 [8.8] 1.4 [4.2] a Calculated using Eqs. 4 – 6 and T = 1250 K [2500 K] b CH3 was excluded from the averages c Calculated for the harmonic reference potential corrected via Eq. 8 d Calculated for the two-harmonic reference potential


Page 27 of 50 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


27 Table 2. Critical point energies (cm–1)a

CH2

1

Inversion Barrier 2069.3

NH3

2

2237.7

CH2CH

2

1826.0

CH3CH3

3

980.1

CH3CHO

3

403.2

CH3OH

3

383.9

CH3OO

3

335.7

CH3CO

3

160.7

3

3

System 3

CH3CH

H2 O2

"b

2660.0

2

Torsional Barrier

Secondary Minimum

115.4 364.6 2593.3

NH2NH2

8

2247.3

747.8 2836.3

CH2OH

4

139.2

1670.3

CH3CH2

6

31.8c

31.8c

CH2CHOH

1

NH2OH

1

HC(O)OH

1

a

538.7

2405.8

1512.2

4345.3

1565.1

CCSD(T)/cc-pVTZ energies without zero point corrections

b c

4632.6

1724.3

Number of classically distinguishable lowest-energy minima

This barrier has mixed torsion/inversion character.



Page 28 of 50

28 References 1

Jasper, A. W.; Hansen, N. Hydrogen-Assisted Isomerizations of Fulvene to Benzene and of Larger Cyclic Aromatic Hydrocarbons. Proc. Combust. Inst. 2013, 34, 279– 287.

2

Klippenstein, S. J. From Theoretical Reaction Dynamics to Chemical Modeling of Combustion. Proc. Comb. Inst. 2017, 36, 77–111.

3

Karton, A.; Ruscic, B.; Martin J. M. L. Benchmark Atomization Energy of Ethane: Importance of Accurate Zero-Point Vibrational Energies and Diagonal BornOppenheimer Corrections for a ‘Simple’ Organic Molecule. THEOCHEM 2007, 811, 345–353.

4

Pfeiffer, F.; Rauhut, G.; Feller, D.; Peterson, K. A. Anharmonic Zero Point Vibrational Energies: Tipping the Scales in Accurate Thermochemistry Calculations. J. Chem. Phys. 2013, 138, 044331.

5

Harding, L. B.; Georgievski, Y.; Klippenstein, S. J. Accurate Anharmonic Zero Point Energies for some Combustion Related Species from Diffusion Monte Carlo, J. Phys. Chem. A, 2017, 121, 4334–4340.

6

McCoy, A. B. Diffusion Monte Carlo Approaches to Investigating the Structure and Vibrational Spectra of Fluxional Systems. Int. Rev. Phys. Chem. 2006, 25, 77–107.

7

McCoy, A. B.; Braams, B. J.; Brown, A.; Huang, X.; Jin, Z.; Bowman, J. M. Ab Initio Diffusion Monte Carlo Calculations of the Quantum Behavior of CH5+ in Full Dimensionality. J. Phys. Chem. A 2004, 23, 4991–4994.

8

McCoy, A. B.; Huang, X.; Carter, S.; Landeweer, M. Y.; Bowman, J. M. Full Dimensional Vibrational Calculations for H5O2+ Using an Ab Initio Potential Energy Surface. J. Chem. Phys. 2005, 122, 061101.

9

Huang, X.; McCoy, A. B.; Bowman, J. M.; Johnson, L. M.; Savage, C.; Dong, F.; Nesbitt, D. J. Quantum Deconstruction of the Infrared Spectrum of CH5+. Science 2006, 311, 60–63.

10

Nielsen, H. H. The Vibration-Rotation Energies of Molecules. Rev. Mod. Phys. 1951, 22, 90–136.

11

Carney, G. D.; Sprandel, L. I.; Kern, C. W. Variational Approaches to the Vibration-


Page 29 of 50 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


29

Rotation Spectroscopy for Polyatomic Molecules. Adv. Chem. Phys. 1978, 37, 305– 379. 12

Bowman, J. M. Self-Consistent Field Energies and Wavefunctions for Coupled Oscillators. J. Chem. Phys. 1978, 68, 608&610.

13

Bowman, J. M.; Christoffel, K.; Tobin, F. Application of SCF-SI Theory to Vibrational Motion in Polyatomic Molecules. J. Phys. Chem. 1979, 83, 905&912.

14

Carter, S.; Bowman, J. M.; Handy, N. C. Extensions and Tests of “Multimode”: A Code to Obtain Accurate Vibration/Rotation Energies of Many-Mode Molecules. Theor. Chem. Acc. 1998, 100, 191–198.

15

Tew, D. P.; Handy, N. C.; Carter, S. A Reaction Surface Hamiltonian Study of Malonaldehyde. J. Chem. Phys. 2006, 125, 084313.

16

Chakraborty, A.; Truhlar, D. G. Converged Vibrational Energy Levels and Quantum Mechanical Vibrational Partition Function of Ethane, J. Chem. Phys. 2006, 124, 184310.

17

Begue, D.; Gohaud, N.; Pouchan, C. A Comparison of Two Methods for Selecting Vibrational Configuration Interaction Spaces on a Heptatomic System: Ethylene Oxide. J. Chem. Phys. 2007, 127, 164115.

18

Keçeli, M.; Shiozaki, T.; Yagi, K.; Hirata, S. Anharmonic Vibrational Frequencies and Vibrationally Averaged Structures of Key Species in Hydrocarbon Combustion: HCO+, HCO, HNO, HOO, HOO–, CH3+, and CH3. Mol. Phys. 2009, 107, 1283–1301.

19

Keçeli, M.; Hirata, S. Size-Extensive Vibrational Self-Consistent Field Method. J. Chem. Phys. 2011, 135, 134108.

20

Yagi, K.; Keçeli, M.; Hirata, S. Optimized Coordinates for Anharmonic Vibrational Structure Theories. J. Chem. Phys. 2012, 137, 204118.

21

Faucheaux, J. A.; Hirata, S. Higher-Order Diagrammatic Vibrational Coupled-Cluster Theory. J. Chem. Phys. 2015, 143, 134105.

22

Ceperley, D. Path Integrals in the Theory of Condensed Helium. Rev. Mod. Phys. 1995, 67, 279–355.

23

Berne, B. J.; Thirumalai, D. On the Simulation of Quantum Systems: Path Integral Methods. Annu. Rev. Phys. Chem. 1986, 37, 401–424.



Page 30 of 50

30

24

Mielke, S. L.; Srinivasan, J.; Truhlar, D. G. Extrapolation and Perturbation Schemes for Accelerating the Convergence of Quantum Mechanical Free Energy Calculations via the Fourier Path-Integral Monte Carlo Method. J. Chem. Phys. 2000, 112, 8758– 8764.

25

Mielke, S. L.;. Truhlar, D. G Displaced-Points Path Integral Method for Including Quantum Effects in the Monte Carlo Evaluation of Free Energies. J. Chem. Phys. 2001, 115, 652–662.

26

Miller, T. F. III; Clary, D. C. Torsional Path Integral Monte Carlo Method for Calculating the Absolute Quantum Free Energy of Large Molecules, J. Chem. Phys. 2003, 119, 68–76.

27

Mielke, S. L.; Truhlar, D. G. Improved Methods for Feynman Path Integral Calculations of Vibrational-Rotational Free Energies and Application to Isotopic Fractionation of Hydrated Chloride Ions. J. Phys. Chem. A 2009, 113, 4817–4827.

28

Mielke, S. L.; Truhlar, D. G. Improved Methods for Feynman Path Integral Calculations and Their Application to Calculate Converged Vibrational–Rotational Partition Functions, Free energies, Enthalpies, Entropies, and Heat Capacities for Methane. J. Chem. Phys. 2015, 142, 044105.

29

Light, J. C.; Hamilton, I. P.; Lill, J. V. Generalized Discrete Variable Approximation in Quantum Mechanics. J. Chem. Phys. 1985, 82, 1400–1409.

30

Tennyson, J. The Calculation of the Vibration-Rotation Energies of Triatomic Molecules Using Scattering Coordinates, Comput. Phys. Rep. 1986, 4, 1–36.

31

Carter, S.; Handy, N. C. A Variational Method for the Determination of the Vibrational (J = 0) Energy Levels of Acetylene, using a Hamiltonian in Internal Coordinates. Comput. Phys. Commun. 1988, 51, 49–58.

32

Chen, R.; Ma, G.; Guo, H. Six-Dimensional Quantum Calculations of Highly Excited Vibrational Energy Levels of Hydrogen Peroxide and its Deuterated Isotopomers. J. Chem. Phys. 2001, 114, 4763–4774.

33

Fabri, C.; Matyus, E.; Furtenbacher, T.; Nemes, L.; Mihaly, B.; Zoltani, T.; Csaszar, A. Variational Quantum Mechanical and Active Database Approaches to the


Page 31 of 50 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


31

Rotational-Vibrational Spectroscopy of Ketene, H2CCO. J. Chem. Phys. 2011, 135, 094307. 34

Christiansen, O. Selected New Developments in Vibrational Structure Theory: Potential Construction and Vibrational Wave Function Calculations. Phys. Chem. Chem. Phys. 2012, 14, 6672&6687.

35

Wodraszka, R.; Manthe, U. Iterative Diagonalization in the Multiconfigurational Time-Dependent Hartree Approach: Ro-vibrational Eigenstates. J. Phys. Chem. A 2013, 117, 7246–7255.

36

Halverson, T.; Poirier, B. Large Scale Exact Quantum Dynamics: Ten Thousand Quantum States of Acetonitrile. Chem. Phys. Lett. 2015, 624, 37–42.

37

Yu, H.-G. An Exact Variational Method to Calculate Rovibrational Spectra of Polyatomic Molecules with Large Amplitude Motion. J. Chem. Phys. 2016, 145, 084109.

38

Carrington, T. Perspective: Computing (Ro-)Vibrational Spectra of Molecules with More than Four Atoms. J. Chem. Phys. 2017, 146, 120902.

39

Avila, G.; Carrington, T. Computing Vibrational Energy Levels of CH4 with a Smolyak Collocation Method. J. Chem. Phys. 2017, 147, 144102.

40

Troe, J. Simplified Models for Anharmonic Numbers and Densities of Vibrational States. I. Application to NO2 and H3+. Chem. Phys. 1995, 190, 381–392.

41

Schmatz, S. Approximate Calculation of Anharmonic Densities of Vibrational States for Very Large Molecules. Chem. Phys. 2008, 346, 198–211.

42

Clabo, D. A.; Allen, W. D.; Remington, R. B.; Yamaguchi, Y.; Schaefer III, H. F. A Systematic Study of Molecular Vibrational Anharmonicity and Vibration-Rotation Interaction by Self-Consistent-Field Higher-Derivative Methods. Asymmetric Top Molecules. Chem. Phys. 1998, 123, 187–239.

43

Schneider, W.; Thiel, W. Anharmonic Force Fields From Analytic Second Derivatives: Method and Application to Methyl Bromide. Chem. Phys. Lett. 1989, 157, 367–373.

44

Wang, F.; Landau, D. P. Efficient, Multiple-Range Random Walk Algorithm to Calculate the Density of States. Phys. Rev. Lett. 2001, 86, 2050–2053.



Page 32 of 50

32

45

Barone, V. Anharmonic Vibrational Properties by a Fully Automated Second-Order Perturbative Approach. J. Chem. Phys. 2005, 122, 014108.

46

Vazquez, J.; Stanton, J. F. Simple(r) Algebraic Equation for Transition Moments of Fundamental Transitions in Vibrational Second-order Perturbation Theory. Mol. Phys. 2006, 104, 377&388.

47

Vazquez, J.; Stanton, J. F. Treatment of Fermi Resonance Effects on Transition Moments in Vibrational Perturbation Theory. Mol. Phys. 2007, 105, 101&109.

48

Nguyen, T. L.; Barker, J. R. Sums and Densities of Fully Coupled Anharmonic Vibrational States: A Comparison of Three Practical Methods. J. Phys. Chem. A 2010, 114, 3718–3730.

49

Pitzer, K. S.; Gwinn, W. D. Energy Levels and Thermodynamic Functions for Molecules with Internal Rotations I. Rigid Frame with Attached Tops. J. Chem. Phys. 1942, 10, 428–440.

50

Truhlar, D. G. A Simple Approximation for the Vibrational Partition Function of a Hindered Internal Rotation. J. Comput. Chem. 1991, 12, 266–270.

51

Knyazev, V. D.; Dubinsky, I. A.; Slagle, I. R.; Gutman, D. Unimolecular Decomposition of t-C4H9 Radical. J. Phys. Chem. 1994, 98, 5279–5289.

52

Forst, W. Sum and Density of States of Polyatomic Systems with Hindered Rotors. J. Computational Chem. 1996, 17, 954–961.

53

McClurg, R. B.; Flagan, R. C.; Goddard III, W. A. The Hindered Rotor Density-ofStates Interpolation Function, J. Chem. Phys. 1997, 106, 6675–6680.

54

Knyazev, V. D. Density of States of One-Dimensional Hindered Internal Rotors and Separability of Rotational Degrees of Freedom. J. Phys. Chem. A 1998, 102, 3916– 3922.

55

Ayala, P. Y.; Schlegel, H. B. Identification and Treatment of Internal Rotation in Normal Mode Vibrational Analysis. J. Chem. Phys. 1998, 108, 2314–2325.

56

Ellingson, B. A.; Lynch, V. A.; Mielke, S. L.; Truhlar, D. G. Statistical Thermodynamics of Bond Torsional Modes: Tests of Separable, Almost-Separable, and Improved Pitzer–Gwinn Approximations, J. Chem. Phys. 2006, 125, 084305.


Page 33 of 50 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


33

57

Sharma, S.; Raman, S.; Green, W. H. Intramolecular Hydrogen Migration in Alkylperoxy and Hydroperoxyalkylperoxy Radicals: Accurate Treatment of Hindered Rotors. J. Phys. Chem. A, 2010, 114, 5689–5701.

58

Zheng, J.; Yu, T.; Papajak, E.; Alecu, I. M.; Mielke, S. L.; Truhlar, D. G. Practical Methods for Including Torsional Anharmonicity in Thermochemical Calculations on Complex Molecules: The Internal-Coordinate Multi-Structural Approximation. Phys. Chem. Chem. Phys. 2011, 13, 10885–10907.

59

Li, Y.-P.; Bell, A. T.; Head-Gordon, M. Thermodynamics of Anharmonic Systems: Uncoupled Mode Approximations for Molecules. J. Chem. Theory Comput. 2016, 12, 2861–2870.

60

Skouteris, D.; Calderini, D.; Barone, V. Methods for Calculating Partition Functions of Molecules Involving Large Amplitude and/or Anharmonic Motions. J. Chem. Theory Comput. 2016, 12, 1011–1018.

61

Doll, J. D. Anharmonic Corrections in Unimolecular Theory: A Monte Carlo Approach. Chem. Phys. Lett. 1980, 72, 139–142.

62

Barker, J. R. Sums of Quantum States for Nonseparable Degrees of Freedom: Multidimensional Monte Carlo Integration. J. Phys. Chem. 1987, 91, 3849–3854.

63

M. Berblinger, C. Schlier, J. Tennyson, S. Miller, Accurate Specific Molecular State Densities by Phase Space Integration. II. Comparison with Quantum Calculations on H3+ and HD2+. J. Chem. Phys. 1992, 96, 6842–6849.

64

Berblinger, M.; Schlier. C. How accurate is the Rice-Ramsperger-Kassel-Marcus theory? The case of H3+. J. Chem. Phys. 1994, 101, 4750–4758.

65

Ming, L.; Nordholm, S.; Schranz, H. W. On the Estimation of Anharmonic Densities of States of Molecules, Chem. Phys. Lett. 1996, 248, 228–236.

66

Taubmann, G.; Schmatz, S. Numbers and Densities of States and Partition Functions From an Efficient Approach to Phase Space Integration. Phys. Chem. Chem. Phys. 2001, 3, 2296–2305.

67

Kamarchik, E.; Jasper, A. W. Anharmonic State Counts and Partition Functions for Molecules Via Classical Phase Space Integrals in Curvilinear Coordinates. J. Chem. Phys. 2013, 138, 194109.



Page 34 of 50

34

68

Kamarchik, E.; Jasper, A. W. Anharmonic Vibrational Properties from Intrinsic nMode State Densities. J. Phys. Chem. Lett. 2013, 4, 2430–2435.

69

Whitten, G. Z.; Rabinovitch, B. S. Accurate and Facile Approximation for Vibrational Energy-Level Sums. J. Chem. Phys. 1963, 38, 2466–2473.

70

Miller, W. H. Classical Path Approximation for the Boltzmann Density Matrix, J. Chem. Phys. 1971, 55, 3146–3156.

71

Wadi, H.; Pollak, E. Accurate Computation of Quantum Densities of States and RRKM Rate Constants for Large Polyatomic Molecules: The STAIR Method. J. Chem. Phys. 1999, 110, 8246–8253.

72

Prudente, F. V.; Varandas, A. J. C. A Direct Evaluation of the Partition Function and Thermodynamic Data for Water at High Temperatures. J. Phys. Chem. A 2002, 106, 6193–6200.

73

Troe, J.; Ushakov, V. G. Anharmonic Rovibrational Numbers and Densities of States for HO2, H2CO, and H2O2. J. Phys. Chem. A 2009, 113, 3940–3945.

74

Kramer, Z. C; Skodje, R. T. A Semiclassical Adiabatic Calculation of State Densities for Molecules Exhibiting Torsion: Application to Hydrogen Peroxide and Its Isotopomers. Theor. Chem. Acc. 2014, 133, 1530 (1-13).

75

Braams, B. J.; Bowman, J. M. Permutationally Invariant Potential Energy Surfaces in High Dimensionality. Int. Rev. Phys. Chem. 2009, 28, 577–606.

76

Lynch, V. A.; Meilke, S. L.; Truhlar, D. G. Accurate Vibrational-Rotational Partition Functions and Standard-State Free Energy Values for H2O2 from Monte Carlo PathIntegral Calculations. J. Chem. Phys. 2004, 121, 5148–5162.

77

McQuarrie, D. A. Statistical Mechanics, University Science Books: 2000.

78

Medvedev, D. M.; Harding. L. B.; Gray, S. K. Methyl Radical: Ab Initio Global Potential Surface, Vibrational Levels and Partition Function, Mol. Phys. 2006, 104, 73–81.

79

Evenhuis, C.; Nyman, G.; Manthe, U. Quantum Dynamics of the CH3 Fragment: A Curvilinear Coordinate System and Kinetic Energy Operators, J. Chem. Phys. 2007, 127, 144302.


Page 35 of 50 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


35

80

Gang, J.; Pilling, M. J.; Robertson, S. H. Partition Functions and Densities of States for Butane and Pentane. J. Chem. Soc. Faraday Trans. 1996, 92, 3509–3518.

81

Robertson, S. H.; Wardlaw, D. M. The Unimolecular Dissociation of the Isopropyl Radical. Chem. Phys. Lett. 1992, 199, 391–396.

82

East, A. L. L.; Radom, L. Ab initio Statistical Thermodynamical Models for the Computation of Third-Law Entropies. J. Chem. Phys. 1997, 106, 6655–6674.

83

Isaacson, A. D.; Truhlar, D. G. The Accuracy of the Pitzer-Gwinn Method for Partition Functions of Anharmonic Vibrational Modes. J. Chem. Phys. 1981, 75, 4090–4094.



Page 36 of 50

36 Figure Captions Fig. 1. (a) Vibrational and (b) rovibrational partition function anharmonicity corrections for the nonfluxional species H2O, HO2, CH2O, CH4, and CH2CH2. Fig. 2. Convergence of the classical vibrational anharmonic correction factor for CH2O using curvilinear coordinates and with respect to the sampling domain parameter ETP = 5000 (green), 10000 (blue), and 15000 (black) cm–1 and Vtest parameters m = 0 and EX = # (dotted), 1 ETP (dashed), and 2 ETP (solid). The result of sampling in Cartesian normal modes with ETP = 15000 cm–1, m = 0, and EX = 2 ETP is shown as a red solid line. Note that for the calculations in this figure Nsamp = 105, whereas elsewhere Nsamp = 106–109. Fig. 3. (a) Vibrational and (b) rovibrational partition function anharmonicity corrections for species with inversions: CH2, CH3, NH3, and CH2CH. Fig 4. Anharmonicity correction to the vibrational density of states #vib of CH2 for all three vibrations (thin solid line) and for the isolated bending coordinate with the remaining coordinates frozen at their equilibrium geometries (thick solid line). The dashed lines show the effect of removing the unphysical large-amplitude contribution to the harmonic reference density of states. The inset shows results for the rovibrational density of states #. Fig. 5. (a) Vibrational and (b) rovibrational partition function anharmonicity corrections for species with symmetric methyl torsions: CH3CH3, CH3OO, CH3OH, CH3CHO, CH3CO, and CH3CH. Fig. 6. (a) Ratio of the one-dimensional hindered rotor partition function to an effective harmonic partition function. (b) Full-dimensional rovibrational anharmonicity correction factors corrected by fHR. Fig. 7. (a) Fitted anharmonic (black) and harmonic (blue) PES as a function of the C–C– H bending and H–C–C–H torsional coordinates. Thin contours are shown for the fitted PES every 100 cm–1. Thick contours for both PESs are shown every 1000 cm–1. (b) Fitted anharmonic (solid black) and harmonic (solid blue) PES along the bending coordinate and for H–C–C–H = 180o. Kinetically separable effective potentials are shown for 600 K (long dashes), 1250 K (medium dashes), and 2500


Page 37 of 50 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


37 K (short dashes). The dotted line shows the magnitude of the Jacobian with an arbitrary scaling factor. Fig. 8. (a) Vibrational and (b) rovibrational partition function anharmonicity corrections for species with asymmetric torsional barriers and inversions: H2O2, NH2NH2, CH2OH, and CH3CH2. Fig. 9. Two dimension contour plot of the NN internal rotation and nitrogen inversion modes of NH2NH2. The inversion angle is defined such that 90° corresponds to a planar nitrogen atom. Fig. 10. Two dimensional contour plot of the NO internal rotation and nitrogen inversion modes of NH2OH. The inversion angle is defined such that 90° corresponds to a planar nitrogen atom. Fig. 11. (a) Vibrational and (b) rovibrational partition function anharmonicity corrections for species with secondary minima: CH2CHOH, HC(O)OH, and NH2OH. Fig. 12. Vibrational anharmonicity corrections relative to the two-harmonic-oscillator models defined by Eq. 12 (thick lines) and Eq. 13 (thin lines).




Page 38 of 50

Page 39 of 50 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 40 of 50

Page 41 of 50 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 42 of 50

Page 43 of 50 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 44 of 50

Page 45 of 50 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 46 of 50

Page 47 of 50 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 48 of 50

Page 49 of 50 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 50 of 50

Anharmonic Rovibrational Partition Functions for Fluxional Species at

Recommend Documents