Polyrate is a suite of computer programs for the calculation of chemical reaction rates of polyatomic species (including atoms and diatoms as special cases) by variational transition state theory (VTST); conventional transition state theory is also supported. Polyrate can calculate the rate constants for both bimolecular reactions and unimolecular reactions, and it can be applied to reactions in the gas phase, liquid solution phase, or solid state and to reactions at gas-solid interfaces. Polyrate can perform VTST calculations on gas-phase reactions with both tight and loose transition states. For tight transition states it uses the reaction-path (RP) variational transition state theory developed by Garrett and Truhlar, and for loose transition states it uses variable-reaction-coordinate (VRC) variational transition state theory developed by Georgievskii and Klippenstein. The RP methods used for tight transition states are conventional transition state theory, canonical variational transition state theory (CVT), and microcanonical variational transition state theory (mu VT) with multidimensional semiclassical approximations for tunneling and nonclassical reflection. For VRC calculations, rate constants may be calculated for canonical or microcanonical ensembles or energy-and total-angular-momentum resolved microcanonical ensembles. Pressure-dependent rate constants for elementary reactions can be computed using system-specific quantum RRK theory (SS-QRRK) with the information obtained from high-pressure limit VTST calculation as input by using the SS-QRRK utility code. Alternatively, Polyrate 2023 may be interfaced with TUMME 2023 for a master-equation treatment of pressure dependence or to obtain phenomenological rate constants for complex mechanisms. Potential energy surfaces may be analytic functions evaluated by subroutines, or they may be implicit surfaces defined by electronic structure input files or interface subroutines containing energies, gradients, and force constants (Hessians). For the latter, Polyrate can be used in conjunction with various interfaces to electronic structure programs for direct dynamics, and it has routines designed to make such interfacing straightforward. Polyrate supports six options for direct dynamics, namely (i) straight single-level direct dynamics, (ii) zero-order interpolated variational transition state theory (IVTST-0), (iii) first-order interpolated variational transition state theory (IVTST-1), (iv) interpolated variational transition state theory by mapping (IVTST-M), (v) variational transition state theory with interpolated single-point energies (VTST-ISPE), and (vi) variational transition state theory with interpolated optimized corrections (VTST-IOC). Polyrate can handle multistructural and torsional-potential anharmonicity in conjunction with the MSTor program. Polyrate 2023 contains 112 test runs, and 46 of these are for direct dynamics calculations; 85 of the test runs are single-level runs, and 27 are dual-level calculations.

A recently developed method for calculating anharmonic vibrational energy levels at nonstationary points along a reaction path that is based on second-order perturbation theory in curvilinear coordinates is combined with variational transition state theory with semiclassical multidimensional tunneling approximations to calculate thermal rate constants for the title reaction. Two different potential energy surfaces were employed for these calculations, an improved version of the author's surface 5 and the WSLFH surface of Wu et al. [J. Chem. Phys. 113, 3150 (2000)]. We present detailed comparisons of rate constants computed for the two surfaces with and without anharmonicity and with various approximations for incorporating tunneling along the reaction path. The results for this system are quite sensitive to the surface employed, the choice of coordinates (curvilinear versus rectilinear), and the inclusion of anharmonicity. A comparison with experiment provides information on the accuracy of these surfaces.

A method for calculating anharmonic vibrational energy levels in asymmetric top and linear systems that is based on second-order perturbation theory in curvilinear coordinates is extended to the bound generalized normal modes at nonstationary points along a reaction path. Explicit formulas for the anharmonicity coefficients, x(ij), and the constant term, E0, are presented, and the necessary modifications for resonance cases are considered. The method is combined with variational transition state theory with semiclassical multidimensional tunneling approximations to calculate thermal rate constants for the HCN/HNC isomerization reaction. Although the results for this system are not very sensitive to the choice of coordinates, we find that the inclusion of anharmonicity leads to a substantial improvement in the vibrational energy levels. We also present detailed comparisons of rate constants computed with and without anharmonicity, with various approximations for incorporating tunneling along the reaction path, and with a more practical approach to calculating the vibrational partition functions needed for larger systems.

Accurate rovibrational levels for HCN up to 24 349 cm−1 above the bottom of the vibrational well and vibrational levels at the saddle point of the HCN/HNC isomerization reaction up to 32 809 cm−1 above the saddle point have been computed and used to obtain partition functions over the temperature range 200–2500 K. Under the rigid-rotor approximation, the rovibrational partition function for HCN is found to be exactly separable into vibrational and rotational contributions to first order. Two approximate approaches, second-order perturbation theory and simple perturbation theory, the latter of which obviates the need for a direct summation over energy levels, are shown to yield vibrational partition functions for HCN and at the saddle point that agree with the accurate values within 2%. In contrast, the usual harmonic approximation leads to errors of up to 24% in the individual partition functions, resulting in differences of between 3.5% and 6.5% in conventional transition state theory rate constants calculated with harmonic versus anharmonic vibrations.

This paper compares different methods for removing resonance effects from second-order perturbation theory calculations of vibrational energies in a variety of systems containing from two to six modes. Both the recently proposed method of Kuhler et al. and the standard approach of Nielsen yield stable energy levels even very close to resonance, with the latter giving smaller average errors in such cases. In addition, the method of Kuhler et al. is observed to affect the ground-state energy, unlike the standard approach. This generally worsens the accuracy of the vibrational partition function at room temperature, especially for systems close to resonance. (C) 1998 American Institute of Physics.

To model proton transfer in biological systems, we consider a modified H5O2+ system, in which each of the two outermost hydrogens (H*) is assigned a large mass in order to represent a backbone. For the potential energy surface of our model, we add a harmonic function, called the "backbone term", to the potential energy function of Ojamae, Shavitt, and Singer for H5O2+. This backbone term holds the H* atoms apart and, thus, provides various oxygen-oxygen distances and barriers for proton transfer. Variational transition-state theory (CVT) rate constants converge for H* masses greater than 1000 amu. These rate constants decrease exponentially as the backbone-backbone equilibrium distance increases. CVT rate constants also decrease as the backbone-term force constant increases and converge in the limit of a large backbone-term force constant. Tunneling effects are more important at low temperature and for larger values of backbone-term force constants or backbone-backbone equilibrium distances. The motion of the system along the minimum energy path from the saddle point to the product involves the motion of the proton between two relatively fixed oxygens followed by fragment motion and relaxation into the product well.

Recent ab initio information of Kraka and Dunning on the reaction OH+H2→H2O+H is used to construct a potential energy surface in the vicinity of the reaction path. The resultant energy surface reproduces the ab initio reactant and product properties and provides a good fit to the ab initio data in the interaction region. Anharmonic vibrational energy levels involving the bound degrees of freedom orthogonal to the reaction coordinate are obtained using perturbation theory through second order for cubic terms and first order for quartic terms, with resonance effects removed. These energy levels are used in the calculation of transmission coefficients and thermal rate constants over the temperature range from 200 to 2400 K. The results are compared with those obtained from harmonic vibrational energy levels.

This paper presents a way of improving second-order perturbation theory calculations by summing contributions of uncoupled excitations to infinite order. For problems involving molecular vibrations, the new theory is shown to give similar results to conventional second-order perturbation theory when the system treated has no near resonances but also to give accurate and stable results even very close to resonance. The new theory is tested by comparison to converged variational calculations for vibrational energy levels of formaldehyde, formaldehyde-d2, and two two-dimensional model subsystems based on formaldehyde.

POLYRATE is a computer program for the calculation of chemical reaction rates of polyatomic species (and also atoms and diatoms as special cases). Version 1.1 was submitted to the CPC Program Library in 1987, and version 4.0.1 was submitted in 1992. Since that time many new capabilities have been added, old ones have been improved, and the code has been made more portable and user-friendly, resulting in the present improved version 6.5. The methods used are variational or conventional transition state theory and multidimensional semiclassical adiabatic and large-curvature approximations for tunneling and nonclassical reflection. Rate constants may be calculated for canonical or microcanonical ensembles or for specific vibrational states of selected modes with translational, rotational, and other vibrational modes treated thermally. Bimolecular and unimolecular reactions and gas-phase, solid-state, and gas-solid interface reactions are all included. Potential energy surfaces may be global analytic functions or implicit functions defined by interpolation from input energies, gradients, and force constants (Hessian matrices) at selected points on a reaction path. The data needed for the dynamics calculations may also be calculated from a global potential energy surface with more accurate calculations at stationary points. The program calculates reaction paths by the Euler, Euler stabilization, or Page-McIver methods. Variational transition states are optimized from among a one-parameter sequence of generalized transition states orthogonal to the reaction path. Tunneling probabilities are calculated by numerical quadrature, using either the centrifugal-dominant-small-curvature approximation, the large-curvature-version-3 approximation, and/or optimized multidimensional tunneling approximations. In the large-curvature case the tunneling probabilities may be summed over final vibrational states for exoergic reactions or initial vibrational states for endoergic reactions.

Recent ab initio information of Kraka and Dunning on the saddle point region for the reaction OH+H2→H2O+H is used to construct an anharmonic potential energy surface valid near the saddle point. Anharmonic vibrational energy levels involving the bound degrees of freedom orthogonal to the reaction coordinate at the saddle point are obtained using perturbation theory through second order for cubic terms and first order for quartic terms, with resonance effects removed. These energy levels are compared to those obtained from an accurate self-consistent field configuration-interaction method, and are used to calculate thermal vibrational partition functions over the temperature range from 200 to 2400 K.

We present a computer program, MORATE (Molecular Orbital RATE calculations), for direct dynamics calculations of unimolecular and bimolecular rate constants of gas-phase chemical reactions involving atoms, diatoms, or polyatomic species. The potential energies, gradients, and higher derivatives of the potential are calculated whenever needed by semiempirical molecular orbital theory without the intermediary of a global or semiglobal fit. The dynamical methods used are conventional or variational transition state theory and multidimensional semiclassical approximations for tunneling and nonclassical reflection. The computer program is conveniently interfaced package consisting of the POLYRATE program, version 4.5.1, for dynamical rate calculations, and the MOPAC program, version 5.03, for semiempirical electronic structure computations. All semiempirical methods available in MOPAC, in particular MINDO/3, MNDO, AM1, and PM3, can be called on to calculate the potential and gradient. Higher derivatives of the potential are obtained by numerical derivatives of the gradient. Variational transition states are found by a one-dimensional search of generalized-transition-state dividing surfaces perpendicular to the minimum-energy path, and tunneling probabilities are evaluated by numerical quadrature.

The energetics of proton transfer in the (H3CH...CH3)- COMPlex, which is one stage in the gas-phase reaction CH4 + CH3- --> CH3- + CH4, have been investigated with ab initio calculations employing a 4-31G basis. The minimum-energy reaction path was determined by a steepest-descent technique in mass-weighted coordinates from the symmetric saddle point (1.35 kcal/mol above separated CH4 + CH3-) to the reactant/product association well (9.25 kcal/mol below separated CH4 + CH3-). Motion along this path is found to occur in two relatively distinct phases: motion of the proton between two fixed carbons followed by separation of the two hydrocarbon fragments. Rate constants computed by variational transition state theory, including an adiabatic approximation for the vibrational modes transverse to the reaction coordinate, demonstrate the importance of including the contributions from tunneling at energies below the top of the 7.88 kcal/mol vibrationally adiabatic barrier for low to moderate temperatures. In addition, incorporating the effects of reaction-path curvature in the tunneling calculation is found to be important, especially at low temperature. However, due to the dual nature of the reaction path, various model barriers fit-to the saddle point information yield rate constants that are much too large at low temperatures, thus underscoring the importance of computing the reaction path explicitly.

POLYRATE is a computer program for the calculation of chemical reaction rates of polyatomic species (and also atoms and diatoms as special cases). Version 1.1 was submitted to the CPC Program Library in 1987, and since that time we have considerably improved the program and made it more portable, and we have added several new capabilities, resulting in the present improved version 4. The methods used are variational or conventional transition state theory and multidimensional semiclassical adiabatic and large-curvature approximations for tunneling and nonclassical reflection. Rate constants may be calculated for canonical or microcanonical ensembles or for specific vibrational states of selected modes with translational, rotational, and other vibrational modes treated thermally. Bimolecular and unimolecular reactions and gas-phase, solid-state, and gas-solid interface reactions are all included. Potential energy surfaces may be global analytic functions or implicit functions defined by interpolation from input energies, gradients, and force constants (Hessian matrices) at selected points on a reaction path. The program calculates reaction paths by the Euler, Euler stabilization, or Page-McIver methods. Variational transition states are optimized from among a one-parameter sequence of generalized transition states orthogonal to the reaction path. Tunneling probabilities are calculated by numerical quadrature, using either the centrifugal-dominat-small-curvature approximation, the large-curvature-version-3 approximation, and / or methods that were available earlier. In the large curvature case the tunneling probabilities may be summed over final vibrational states for exoergic reactions or initial vibrational states for endoergic reactions.

Two new potential energy surfaces for modeling the dynamics of the reaction OH + H2 --> H2O + H are presented. While the form of these surfaces is relatively simple, varying the surface parameters along the reaction path allows these surfaces to reproduce the ab initio information of Walch and Dunning as well as either the barrier shape computed by Dunning, Harding, and Kraka or the experimental rate constants of Ravishankara et al. By comparing the results obtained from these surfaces and an earlier one by Schatz and Elgersma, we are able to study the sensitivity of the predicted dynamics on both the form of the surface and the data to which the surface is fit. In particular, the results indicate that both the barrier shape and the degree of reaction path curvature can strongly influence the predicted rate constants at lower temperatures.

Positions and widths for the lowest 1-SIGMA-g+ doubly excited autoionizing states of H-2 at several internuclear separations have been obtained by the calculation of Siegert eigenvalues. This approach involves the direct computation of the complex resonance energy in a basis set of both real and complex Slater orbitals. When the complex orbitals are centered between the two atoms, the numerical results are in fair agreement with previous theoretical treatments, although the present width does not rise quite as much for larger internuclear separations. In addition, the sensitivity of the present results on the basis set is studied and improvements to the basis set which should provide converged positions and widths for these resonances are discussed.

We propose and test a very simple method for calculating equilibrium constants from quartic force fields.

Positions and widths for the lowest 1Σ+g doubly excited autoionizing states of H2 at several internuclear separations have been obtained by the calculation of Siegert eigenvalues. This approach involves the direct computation of the complex resonance energy in a basis set of both real and complex Slater orbitals. When the complex orbitals are centered between the two atoms, the numerical results are in fair agreement with previous theoretical treatments, although the present width does not rise quite as much for larger internuclear separations. In addition, the sensitivity of the present results on the basis set is studied and improvements to the basis set which should provide converged positions and widths for these resonances are discussed.

Several variational methods are applied to the calculation of the position and width of the lowest 1S resonance state of H−, which is the simplest physical example of an electronic Feshbach resonance. These methods include two different versions of the analytic continuation of stabilization graphs that enforce the correct branch-point structure and two versions of the complex-stabilization approach, one that stabilizes the complex resonance energy with respect to the exponents of the complex orbital(s) and one that stabilizes it with respect to both the real and complex orbital exponents. The calculations involve medium-, large-, and very-large-sized basis sets of Gaussian orbitals and full configuration interaction (CI). The use of the same basis sets with the various methods allows for detailed comparisons among them. Although the sensitivity of the results to the fit parameters prevents true convergence, reliable estimates of the position and width of this resonance (about four-figure accuracy in the position and two-figure accuracy in the width) are obtained both from a version of the analytic continuation of stabilization graphs that employs one eigenvalue of a real, Hermitian Hamiltonian matrix but enforces the correct branch-point structure and from a complex-stabilization approach that involves complex basis functions and a non-Hermitian Hamiltonian matrix. In the former approach, we find that the results are less accurate when two eigenvalues of the Hamiltonian matrix are employed in the analytic continuation, possibly due to interactions with excited resonance states. For the latter approach, we show that good results can be obtained with basis sets containing a single complex orbital if the resonance energy is also stabilized with respect to an analytic continuation of the real orbital exponents, but that there is no advantage in using two complex orbitals with close exponents.

Purely vibrational energy levels and partition functions are calculated using three different potential energy surfaces for the H2O molecule. Results obtained with perturbation-theory, independent-normal-mode (INM), and harmonic approximations are compared with accurate values. For the cases considered here, the expected improvement that perturbation theory provides over the corresponding harmonic treatment is found to be substantial, while the INM approximation leads to results which are worse than the corresponding harmonic ones. In fact, we show that reliable partition functions for these potential surfaces can be obtained when resonance contributions are removed from the perturbation-theory treatment, and we propose a theoretical criterion for deciding when a particular interaction should be treated as resonant.