The cost of calculating second derivatives of the energy, or nuclear hessians, in the course of quantum chemical analyses can be prohibitive for systems containing hundreds of atoms. In particular, when searching for reaction transition states (TSs), only a few eigenvalues and eigenvectors, and not the full hessian, are required. Here, a method is described that can eliminate the need for hessian calculations for both TS searches as well as characterization of stationary points. A finite differences implementation of the Davidson method that uses only first derivatives of the energy to calculate the lowest eigenvalues and eigenvectors of the hessian is discussed. When implemented in conjunction with a double-ended interpolation method for generating TS guesses, such as the freezing string method (FSM), an approximate hessian can be constructed in lieu of the full hessian as input to any quasi-Newton TS optimization routine. With equal ease, the finite differences Davidson approach can be implemented at the end of geometry optimization for verifying stationary points on a potential energy surface. The approach scales one power of system size lower than exact hessian calculation since the rate of convergence is approximately independent of the size of the system. Therefore, it achieves significant cost savings relative to exact hessian calculation when applied to both stationary point characterization as well as TS search, particularly when analytical hessians are not available or require substantial computational effort.
The TS search approach is a useful tool for reaction kinetics and catalysis studies. Zeolite catalysts are employed extensively in industry owing to their high Brønsted acidity and shape selective properties, which are probed typically using monomolecular cracking and dehydrogenation reactions of alkanes. The TS search method is combined with hybrid quantum mechanics/molecular mechanics (QM/MM), and a modified harmonic oscillator approximation in order to calculate intrinsic activation parameters for monomolecular reactions of n-butane. The first study calculates TSs for all cracking and dehydrogenation pathways in MFI. Based on an examination of adsorption enthalpies and intrinsic activation energies for these reactions at active sites located at the channel intersection as well as the sinusoidal channel in MFI, the analysis concludes that reaction energetics are highly sensitive to the active site location due to varying acidities and non-bonding framework-substrate interactions.
The second investigation extends the QM/MM approach to examine the sensitivity of intrinsic reaction kinetics to zeolite pore topology. Monomolecular cracking and dehydrogenation reactions of n-butane are examined in six zeolite frameworks - TON, SVR, MFI, MEL, STF and MWW, with active sites located within channels, channel intersections and cage geometries. By analyzing calculated intrinsic enthalpies and entropies of activation together with experimental values, the sensitivity of cracking and dehydrogenation pathways to active site location is examined for all site types. Dehydrogenation exhibits a surprising preference for the methyl pathway in cages in spite of the higher barrier relative to methylene, which points towards significant entropy compensation occurring at these active sites. However, although computed enthalpies of activation are in good agreement with experiment, thermochemical approximations that better account for anharmonic contributions are required to accurately determine entropy differences between these pathways.
The hessian-free finite differences Davidson approach can also be extended to the space of molecular orbital coefficients. Wavefunction stability analysis is commonly applied to converged self-consistent field (SCF) solutions to verify whether the electronic energy is a local minimum with respect to second-order variation in the orbitals, by calculating the lowest eigenvalue of the electronic hessian. Analytical expressions for the electronic hessian are unavailable for some advanced post-Hartree–Fock (HF) wave function methods and even certain Kohn–Sham (KS) density functionals. Calculating full finite difference hessians for even small molecules can prove intractable in such cases. To address this issue, the hessian-vector product within the Davidson scheme is formulated as a finite difference of the electronic gradient with respect to orbital perturbations. As a model application, following the lowest eigenvalue of the orbital-optimized second-order Møller–Plesset perturbation theory (OOMP2) hessian during H2 dissociation reveals the surprising stability of the spin-restricted solution at all separations, with a second independent unrestricted solution. A single stable solution can be recovered by using the regularized OOMP2 method (δ-OOMP2), which contains a level shift. Internal and external stability analyses are also performed for SCF solutions of a recently developed range-separated hybrid density functional, ωB97X-V, for which the analytical hessian is not yet available due to the complexity of its long-range non-local VV10 correlation functional.
