Interpolating moving least-squares methods for fitting potential energy surfaces: a strategy for efficient automatic data point placement in high dimensions

An accurate and efficient method for automated molecular global potential energy surface (PES) construction and fitting is demonstrated. An interpolating moving least-squares (IMLS) method is developed with the flexibility to fit various ab initio data: (1) energies, (2) energies and gradients, or (3) energies, gradients, and Hessian data. The method is automated and flexible so that a PES can be optimally generated for trajectories, spectroscopy, or other applications. High efficiency is achieved by employing local IMLS in which fitting coefficients are stored at a limited number of expansion points, thus eliminating the need to perform weighted least-squares fits each time the potential is evaluated. An automatic point selection scheme based on the difference in two successive orders of IMLS fits is used to determine where new ab initio data need to be calculated for the most efficient fitting of the PES. A simple scan of the coordinate is shown to work well to identify these maxima in one dimension, but this search strategy scales poorly with dimension. We demonstrate the efficacy of using conjugate gradient minimizations on the difference surface to locate optimal data point placement in high dimensions. Results that are indicative of the accuracy, efficiency, and scalability are presented for a one-dimensional model potential (Morse) as well as for three-dimensional (HCN), six-dimensional (HOOH), and nine-dimensional (CH4) molecular PESs.     

Interpolating moving least-squares methods for fitting potential energy surfaces: an application to the H2CN unimolecular reaction

Classical trajectories have been used to compute rates for the unimolecular reaction H2CN-->H+HCN on a fitted ab initio potential energy surface (PES). The ab initio energies were obtained from CCSD(T)/aug-cc-pvtz electronic structure calculations. The ab initio energies were fitted by the interpolating moving least-squares (IMLS) method. This work continues the development of the IMLS method for producing ab initio PESs for use in molecular dynamics simulations of many-atom systems. A dual-level scheme was used in which the preliminary selection of data points was done using a low-level theory and the points used for fitting the final PES were obtained at the desired higher level of theory. Classical trajectories were used on various low-level IMLS fits to tune the fit to the unimolecular reaction under study. Procedures for efficiently picking data points, selecting basis functions, and defining cutoff limits to exclude distant points were investigated. The accuracy of the fitted PES was assessed by comparing interpolated values of quantities to the corresponding ab initio values. With as little as 330 ab initio points classical trajectory rate constants were converged to 5%-10% and the rms error over the six-dimensional region sampled by the trajectories was a few tenths of a kcal/mol.     

Interpolating moving least-squares methods for fitting potential-energy surfaces: further improvement of efficiency via cutoff strategies

In standard applications of interpolating moving least squares (IMLS) for fitting a potential-energy surface (PES), all available ab initio points are used. Because remote ab initio points negligibly influence IMLS accuracy and increase IMLS time-to-solution, we present two methods to locally restrict the number of points included in a particular fit. The fixed radius cutoff (FRC) method includes ab initio points within a hypersphere of fixed radius. The density adaptive cutoff (DAC) method includes points within a hypersphere of variable radius depending on the point density. We test these methods by fitting a six-dimensional analytical PES for hydrogen peroxide. Both methods reduce the IMLS time-to-solution by about an order of magnitude relative to that when no cutoff method is used. The DAC method is more robust and efficient than the FRC method.     

Interpolating moving least-squares methods for fitting potential energy surfaces: Analysis of an application to a six-dimensional system

The basic formal and numerical aspects of different degree interpolated moving least-squares (IMLS) methods are applied to a six-dimensional potential energy surface (PES) of the HOOH molecule, for which an analytic ("exact") potential is available in the literature. The results of systematic investigations of the effects of weight function parameters, the degree and partial degree of IMLS, the number of data points allowed, and the optimal automatic point selection of data points up to full third-degree IMLS fits are reported. With partial reduction of cross terms and automatic point selection the full six-dimensional HOOH PES can be fit over a range of 100 kcal/mol to an accuracy of less than 1 kcal/mol with approximately 1350 ab initio points.     

Potential energy surface of HDO up to 25,000 cm-1

A new spectroscopically determined potential energy surface (PES) for HD(16)O is presented. This surface is constructed by adjusting the high accuracy ab initio PES of Polyansky et al. [Science 299, 539 (2003)] by fitting to both published experimental data and our still unpublished data. This refinement used experimentally derived term values up to 25,000 cm(-1) and with J< or =8: a data set of 3478 energy levels once some levels with ambiguous assignment is excluded. To improve the extrapolation properties of the empirical PES, the restraint that the resulting PESs remain close to the ab initio surface was imposed. The new HDO_07 PES reproduces the experimental data, including high J levels not included in the fit, with a root mean square error of 0.035 cm(-1). Predictions for rotation-vibration term values up to J=12 are made.     

Gradient incorporation in one-dimensional applications of interpolating moving least-squares methods for fitting potential energy surfaces

We present several approaches to use gradients in higher degree interpolating moving least squares (IMLS) methods for representing a potential energy surface (PES). General procedures are developed to obtain smooth approximations of the PES and its derivatives from quasi-uniform sets of energy and gradient data points. These methods are illustrated and analyzed for the Morse oscillator and a 1-D slice of the ground-state PES for the HCO radical computed using density functional theory. Variations in the IMLS fits with the number and distribution of points and the degree of the polynomial fitting basis set are examined. We determine the effects of gradient inclusion on the accuracy of the IMLS values of the energy, first and second derivatives for two 1-D test cases. Gradient inclusion reduces the number of data points required by up to 40%.     

High-dimensional ab initio potential energy surfaces for reaction dynamics calculations

There has been great progress in the development of potential energy surfaces (PESs) for reaction dynamics that are fits to ab initio energies. The fitting techniques described here explicitly represent the invariance of the PES with respect to all permutations of like atoms. A review of a subset of dynamics calculations using such PESs (currently 16 such PESs exist) is then given. Bimolecular reactions of current interest to the community, namely, H + CH(4) and F + CH(4), are focused on. Unimolecular reactions are then reviewed, with a focus on the photodissociation dynamics of H(2)CO and CH(3)CHO, where so-called "roaming" pathways have been discovered. The challenges for electronically non-adiabatic reactions, and associated PESs, are presented with a focus on the OH* + H(2) reaction. Finally, some thoughts on future directions and challenges are given.     

CVRQD ab initio ground-state adiabatic potential energy surfaces for the water molecule

The high accuracy ab initio adiabatic potential energy surfaces (PESs) of the ground electronic state of the water molecule, determined originally by Polyansky et al. [Science 299, 539 (2003)] and called CVRQD, are extended and carefully characterized and analyzed. The CVRQD potential energy surfaces are obtained from extrapolation to the complete basis set of nearly full configuration interaction valence-only electronic structure computations, augmented by core, relativistic, quantum electrodynamics, and diagonal Born-Oppenheimer corrections. We also report ab initio calculations of several quantities characterizing the CVRQD PESs, including equilibrium and vibrationally averaged (0 K) structures, harmonic and anharmonic force fields, harmonic vibrational frequencies, vibrational fundamentals, and zero-point energies. They can be considered as the best ab initio estimates of these quantities available today. Results of first-principles computations on the rovibrational energy levels of several isotopologues of the water molecule are also presented, based on the CVRQD PESs and the use of variational nuclear motion calculations employing an exact kinetic energy operator given in orthogonal internal coordinates. The variational nuclear motion calculations also include a simplified treatment of nonadiabatic effects. This sophisticated procedure to compute rovibrational energy levels reproduces all the known rovibrational levels of the water isotopologues considered, H(2) (16)O, H(2) (17)O, H(2) (18)O, and D(2) (16)O, to better than 1 cm(-1) on average. Finally, prospects for further improvement of the ground-state adiabatic ab initio PESs of water are discussed.     

Nitrous oxide dimer: a new potential energy surface and rovibrational spectrum of the nonpolar isomer

The spectrum of nitrous oxide dimer was investigated by constructing new potential energy surfaces using coupled-cluster theory and solving the rovibrational Schro?dinger equation with a Lanczos algorithm. Two four-dimensional (rigid monomer) global ab initio potential energy surfaces (PESs) were made using an interpolating moving least-squares (IMLS) fitting procedure specialized to describe the interaction of two linear fragments. The first exploratory fit was made from 1646 CCSD(T)/3ZaP energies. Isomeric minima and connecting transition structures were located on the fitted surface, and the energies of those geometries were benchmarked using complete basis set (CBS) extrapolations, counterpoise (CP) corrections, and explicitly correlated (F12b) methods. At the geometries tested, the explicitly correlated F12b method produced energies in close agreement with the estimated CBS limit. A second fit to 1757 data at the CCSD(T)-F12b/VTZ-F12 level was constructed with an estimated fitting error of less than 1.5?cm(-1). The second surface has a global nonpolar O-in minimum, two T-shaped N-in minima, and two polar minima. Barriers between these minima are small and some wave functions have amplitudes in several wells. Low-lying rovibrational wave functions and energy levels up to about 150?cm(-1) were computed on the F12b PES using a discrete variable representation/finite basis representation method. Calculated rotational constants and intermolecular frequencies are in very close agreement with experiment.     

Interpolating moving least-squares methods for fitting potential energy surfaces: Improving efficiency via local approximants

The local interpolating moving least-squares (IMLS) method for constructing potential energy surfaces is investigated. The method retains the advantageous features of the IMLS approach in that the ab initio derivatives are not required and high degree polynomials can be used to provide accurate fits, while at the same time it is much more efficient than the standard IMLS approach because the least-squares solutions need to be calculated only once at the data points. Issues related to the implementation of the local IMLS method are investigated and the accuracy is assessed using HOOH as a test case. It is shown that the local IMLS method is at the same level of accuracy as the standard IMLS method. In addition, the scaling of the method is found to be a power law as a function of number of data points N, N(-q). The results suggest that when fitting only to the energy values for a d-dimensional system by using a Qth degree polynomial the power law exponent q approximately Qd when the energy range fitted is large (e.g., E<100 kcalmol for HOOH), and q>Qd when the energy range fitted is smaller (E<30 kcalmol) and the density of data points is higher. This study demonstrates that the local IMLS method provides an efficient and accurate means for constructing potential energy surfaces.     

High-level ab initio potential energy surfaces and vibrational energies of H2CS

Six-dimensional (6D) potential energy surfaces (PESs) of H(2)CS have been generated ab initio using the recently proposed explicitly correlated (F12) singles and doubles coupled cluster method including a perturbational estimate of connected triple excitations, CCSD(T)-F12b [T. B. Adler, G. Knizia, and H.-J. Werner, J. Chem. Phys. 127, 221106 (2007)] in conjunction with F12-optimized correlation consistent basis sets. Core-electron correlation, high-order correlation, scalar relativistic, and diagonal Born-Oppenheimer terms were included as additive high-level (HL) corrections. The resulting 6D PESs were represented by analytical functions which were used in variational calculations of the vibrational term values below 5000 cm(-1). The best PESs obtained with and without the HL corrections, VQZ-F12(*HL) and VQZ-F12?, reproduce the fundamental vibrational wavenumbers with mean absolute deviations of 1.13 and 1.22 cm(-1), respectively. A detailed analysis of the effects of the HL corrections shows how the VQZ-F12 results benefit from error cancellation. The present purely ab initio PESs will be useful as starting points for empirical refinements towards an accurate "spectroscopic" PES of H(2)CS.     

Molecular potential energy surfaces constructed from interpolation of systematic fragment surfaces

A systematic method for approximating the ab initio electronic energy of molecules from the energies of molecular fragments has previously been presented. Here it is shown that this approach provides a feasible, systematic method for constructing a global molecular potential energy surface (PES) for reactions of a moderate-sized molecule from the corresponding surfaces for small molecular fragments. The method is demonstrated by construction of PESs for the reactions of a hydrogen atom with propane and n-pentane.     

Application of the modified Shepard interpolation method to the determination of the potential energy surface for a molecule-surface reaction: H2 + Pt(111)

We have used a modified Shepard (MS) interpolation method, initially developed for gas phase reactions, to build a potential energy surface (PES) for studying the dissociative chemisorption of H2 on Pt(111). The aim was to study the efficiency and the accuracy of this interpolation method for an activated multidimensional molecule-surface reactive problem. The strategy used is based on previous applications of the MS method to gas phase reactions, but modified to take into account special features of molecule-surface reactions, like the presence of many similar reaction pathways which vary only slightly with surface site. The efficiency of the interpolation method was tested by using an already existing PES to provide the input data required for the construction of the new PES. The construction of the new PES required half as many ab initio data points as the construction of the old PES, and the comparison of the two PESs shows that the method is able to reproduce with good accuracy the most important features of the H2 + Pt(111) interaction potential. Finally, accuracy tests were done by comparing the results of dynamics simulations using the two different PESs. The good agreement obtained for reaction probabilities and probabilities for rotationally and diffractionally inelastic scattering shows clearly that the MS interpolation method can be used efficiently to yield accurate PESs for activated molecule-surface reactions.     

A new "spectroscopic" potential energy surface for formaldehyde in its ground electronic state

We report a new "spectroscopic" potential energy surface (PES) of formaldehyde (H(2)(12)C(16)O) in its ground electronic state, obtained by refining an ab initio PES in a least-squares fitting to the experimental spectroscopic data for formaldehyde currently available in the literature. The ab initio PES was computed using the CCSD(T)/aug-cc-pVQZ method at 30 840 geometries that cover the energy range up to 44 000 cm(-1) above equilibrium. Ro-vibrational energies of formaldehyde were determined variationally for this ab initio PES by means of the program TROVE [Theoretical ROtation-Vibration Energies; S. N. Yurchenko, W. Thiel, and P. Jensen, J. Mol. Spectrosc. 245, 126 (2007)]. The parameter values in the analytical representation of the PES were optimized in fittings to 319 ro-vibrational energies with J = 0, 1, 2, and 5. The initial parameter values in the fittings were those of the ab initio PES, the ro-vibrational eigenfunctions obtained from this PES served as a basis set during the fitting process, and constraints were imposed to ensure that the refined PES does not deviate unphysically from the ab initio one in regions of configuration space not sampled by the experimental data. The resulting refined PES, referred to as H(2)CO-2011, reproduces the available experimental J ≤ 5 data with a root-mean-square error of 0.04 cm(-1).     

Interpolating moving least-squares methods for fitting potential energy surfaces: applications to classical dynamics calculations

As a continuation of our efforts to develop efficient and accurate interpolating moving least-squares (IMLS) methods for generating potential energy surfaces, we carry out classical trajectories and compute kinetics properties on higher degree IMLS surfaces. In this study, we have investigated the choice of coordinate system, the range of points (i.e., the cutoff radius) used in fitting, and strategies for selections of data points and basis elements. We illustrate and test the method by applying it to hydrogen peroxide (HOOH). In particular, reaction rates for the O-O bond breaking in HOOH are calculated on fitted surfaces using the classical trajectory approach to test the accuracy of the IMLS method for providing potentials for dynamics calculations.     

Automatic generation of potential energy and property surfaces of polyatomic molecules in normal coordinates

A procedure for the automatic construction of Born-Oppenheimer (BO) potential energy and molecular property surfaces in rectilinear normal coordinates is presented and its suitability and accuracy when combined with vibrational structure calculations are assessed. The procedure relies on a hierarchical n-mode representation of the BO potential energy or molecular property surface, where the n-mode term of the sequence of potentials/molecular properties includes only the couplings between n or less vibrational degrees of freedom. Each n-mode cut of the energy/molecular property surface is first evaluated in a grid of points with ab initio electronic structure methods. The ab initio data are then spline interpolated and a subsequent polynomial fitting provides an analytical semiglobal representation for use in vibrational structure programs. The implementation of the procedure is outlined and the accuracy of the method is tested on water and difluoromethane. Strategies for improving the proposed algorithm are also discussed.     

Calibration-quality adiabatic potential energy surfaces for H3(+) and its isotopologues

Calibration-quality ab initio adiabatic potential energy surfaces (PES) have been determined for all isotopologues of the molecular ion H(3)(+). The underlying Born-Oppenheimer electronic structure computations used optimized explicitly correlated shifted Gaussian functions. The surfaces include diagonal Born-Oppenheimer corrections computed from the accurate electronic wave functions. A fit to the 41,655 ab initio points is presented which gives a standard deviation better than 0.1 cm(-1) when restricted to the points up to 6000 cm(-1) above the first dissociation asymptote. Nuclear motion calculations utilizing this PES, called GLH3P, and an exact kinetic energy operator given in orthogonal internal coordinates are presented. The ro-vibrational transition frequencies for H(3)(+), H(2)D(+), and HD(2)(+) are compared with high resolution measurements. The most sophisticated and complete procedure employed to compute ro-vibrational energy levels, which makes explicit allowance for the inclusion of non-adiabatic effects, reproduces all the known ro-vibrational levels of the H(3)(+) isotopologues considered to better than 0.2 cm(-1). This represents a significant (order-of-magnitude) improvement compared to previous studies of transitions in the visible. Careful treatment of linear geometries is important for high frequency transitions and leads to new assignments for some of the previously observed lines. Prospects for further investigations of non-adiabatic effects in the H(3)(+) isotopologues are discussed. In short, the paper presents (a) an extremely accurate global potential energy surface of H(3)(+) resulting from high accuracy ab initio computations and global fit, (b) very accurate nuclear motion calculations of all available experimental line data up to 16,000 cm(-1), and (c) results suggest that we can predict accurately the lines of H(3)(+) towards dissociation and thus facilitate their experimental observation.     

Global potential energy surface, vibrational spectrum, and reaction dynamics of the first excited (Ã (2)A(')) state of HO(2)

The authors report extensive high-level ab initio studies of the first excited (A??(2)A(')) state of HO(2). A global potential energy surface (PES) was developed by spline-fitting 17?000 ab initio points at the internal contracted multireference configuration interaction (icMRCI) level with the AVQZ basis set. To ascertain the spectroscopic accuracy of the PES, the near-equilibrium region of the molecule was also investigated using three interpolating moving least-squares-based PESs employing dynamically weighted icMRCI methods in the complete basis set limit. Vibrational energy levels on all four surfaces agree well with each other and a new assignment of some vibrational features is proposed. In addition, the dynamics of both the forward and reverse directions of the H+O(2)(a??(1)Δ(g))?OH+O reaction (J=0) were studied using an exact wave packet method. The reactions are found to be dominated by sharp resonances.     

Spin-orbit coupled potential energy surfaces and properties using effective relativistic coupling by asymptotic representation

A new method has been reported recently [H. Ndome, R. Welsch, and W. Eisfeld, J. Chem. Phys. 136, 034103 (2012)] that allows the efficient generation of fully coupled potential energy surfaces (PESs) including derivative and spin-orbit (SO) coupling. The method is based on the diabatic asymptotic representation of the molecular fine structure states and an effective relativistic coupling operator and therefore is called effective relativistic coupling by asymptotic representation (ERCAR). The resulting diabatic spin-orbit coupling matrix is constant and the geometry dependence of the coupling between the eigenstates is accounted for by the diabatization. This approach allows to generate an analytical model for the fully coupled PESs without performing any ab initio SO calculations (except perhaps for the atoms) and thus is very efficient. In the present work, we study the performance of this new method for the example of hydrogen iodide as a well-established test case. Details of the diabatization and the accuracy of the results are investigated in comparison to reference ab initio calculations. The energies of the adiabatic fine structure states are reproduced in excellent agreement with reference ab initio data. It is shown that the accuracy of the ERCAR approach mainly depends on the quality of the underlying ab initio data. This is also the case for dissociation and vibrational level energies, which are influenced by the SO coupling. A method is presented how one-electron operators and the corresponding properties can be evaluated in the framework of the ERCAR approach. This allows the computation of dipole and transition moments of the fine structure states in good agreement with ab initio data. The new method is shown to be very promising for the construction of fully coupled PESs for more complex polyatomic systems to be used in quantum dynamics studies.     

Exploring the accuracy level of new potential energy surfaces for the F + HD reactions: from exact quantum rate constants to the state-to-state reaction dynamics

Exact quantum reactive scattering calculations in the collision energy range 1-250 meV have been carried out for both the isotopic product channels of the title system. The dynamical studies compares an ab initio potential energy surface (PES) recently appeared in the literature (J. Chem. Phys., 2008, 129, 011103) with other phenomenological PESs. Vibrational branching ratios, cross sections and rate constants are presented and compared with molecular beam scattering experiments as well as with chemical kinetics data. In particular, the agreement of the vibrational branching ratios with experimental measurements is improved with respect to previous studies on other PESs, mainly because of the presence of a broad peak in the HF(v' = 3) integral cross section completely absent in the previous simulations. This feature, observed by molecular beam experiments, is the fingerprint of a new reaction mechanism operative in the dynamics described by the new PES. A conjecture for its origin, able to explain many of its characteristic aspects, is analyzed and discussed.