In reaction path optimization, such as the calculation of a minimum energy path (MEP) between given reactant and product configurations of atoms, it is advantageous to start with an initial guess where the close proximity of atoms is avoided and bonds are not unnecessarily broken only to be reformed later. When the configurations of the atoms are described with Cartesian coordinates, a linear interpolation between the end points can be problematic, and a better option is provided by the so-called image dependent pair potential (IDPP) approach where interpolated pairwise distances are generated to form an objective function that can be used to construct an improved initial path. When started with a linear interpolation, this method can, however, still lead to unnecessary bond breaking in, for example, reactions in which a molecular subgroup undergoes significant rotation. In the method presented here, this problem is addressed by constructing the path gradually, introducing images sequentially starting from the vicinity of the end points while the distance between images in the central region is larger. The distribution of images is controlled by systematically scaling the tightness of springs acting between the images until the desired number of images is obtained, and they are evenly spaced. This procedure generates an initial path on the IDPP surface, a task that requires negligible computational effort, as no evaluation of the energy of the system is needed. The calculation of the MEP, typically using electronic structure calculations, is then subsequently carried out in a way that makes efficient use of parallel computing with the nudged elastic band method. Several examples of reactions are given where the linear interpolation IDPP (LI-IDPP) method yields problematic paths with unnecessary bond breaking in some of the intermediate images, while the sequential IDPP (S-IDPP) method yields paths that are significantly closer to realistic MEPs.

Recent Rydberg spectroscopy measurements of a diamine molecule, N,N′-dimethylpiperazine (DMP), indicate the existence of a localized electronic state as well as a delocalized electronic state. This implies that the cation, DMP+, can similarly have its positive charge either localized on one of the N atoms or delocalized over both. This interpretation of the experiments has, however, been questioned based on coupled cluster calculations. In this article, results of high-level multireference configuration interaction calculations are presented where a localized state of DMP+ is indeed found to be present with an energy barrier separating it from the delocalized state. The energy difference between the two states is in excellent agreement with the experimental estimate. The results presented here, therefore, support the original interpretation of the experiments and illustrate a rare shortcoming of CCSD(T), the “gold standard” of quantum chemistry. These results have implications for the development of density functionals, as most functionals fail to produce the localized state.

The climbing image nudged elastic band method (CI-NEB) is used to identify reaction coordinates and to find saddle points representing transition states of reactions. It can make efficient use of parallel computing as the calculations of the discretization points, the so-called images, can be carried out simultaneously. In typical implementations, the images are distributed evenly along the path by connecting adjacent images with equally stiff springs. However, for systems with a high degree of flexibility, this can lead to poor resolution near the saddle point. By making the spring constants increase with energy, the resolution near the saddle point is improved. To assess the performance of this energy-weighted CI-NEB method, calculations are carried out for a benchmark set of 121 molecular reactions. The performance of the method is analyzed with respect to the input parameters. Energy-weighted springs are found to greatly improve performance and result in successful location of the saddle points in less than a thousand energy and force evaluations on average (about a hundred per image) using the same set of parameter values for all of the reactions. Even better performance is obtained by stopping the calculation before full convergence and complete the saddle point search using an eigenvector following method starting from the location of the climbing image. This combination of methods, referred to as NEB-TS, turns out to be robust and highly efficient as it reduces the average number of energy and force evaluations down to a third, to 305. An efficient and flexible implementation of these methods has been made available in the ORCA software.

The relative stability of a localized and delocalized electronic state in the same molecule, the N,N' -Dimethylpiperazine cation, is calculated at various levels of theory up to multireference configuration interaction (MRCI+Q). This system has received a great deal of attention because of recent experimental studies of corresponding Rydberg states of the molecule and the failure of most density functional approximations to produce a metastable localized state. A cut through the energy surface involving two dihedral angles is generated at the level of MRCI+Q as well as Hartree-Fock (HF), Möller-Plesset second order perturbation theory (MP2), coupled cluster theory with and without perturbative triple excitations (CCSD and CCSD(T)) and complete active space self-consistent field calculations with and without perturbative correction (CASSCF and NEVPT2). Remarkably, while CCSD produces a localized state, CCSD(T) does not, and similarly, large active-space CASSCF does while NEVPT2 does not. The inclusion of dynamic correlation in a perturbative way thus adversely affects the accuracy of the calculation, most notably for CCSD(T), the 'golden standard'. The MRCI+Q results are in close correspondence with the experimental results, as well as CAS(19,20) DMRG-CASSCF calculations. The results presented here establish a benchmark system for the study of electronic state localization.

The minimum mode following method can be used to find saddle points on an energy surface by following a direction guided by the lowest curvature mode. Such calculations are often started close to a minimum on the energy surface to find out which transitions can occur from an initial state of the system, but it is also common to start from the vicinity of a first order saddle point making use of an initial guess based on intuition or more approximate calculations. In systems where accurate evaluations of the energy and its gradient are computationally intensive, it is important to exploit the information of the previous evaluations to enhance the performance. Here, we show that the number of evaluations required for convergence to the saddle point can be significantly reduced by making use of an approximate energy surface obtained by a Gaussian process model based on inverse inter-atomic distances, evaluating accurate energy and gradient at the saddle point of the approximate surface and then correcting the model based on the new information. The performance of the method is tested with start points chosen randomly in the vicinity of saddle points for dissociative adsorption of an H2 molecule on the Cu(110) Surface and three gas phase chemical reactions.

Calculations of minimum energy paths for atomic rearrangements using the nudged elastic band method can be accelerated with Gaussian process regression to reduce the number of energy and atomic force evaluations needed for convergence. Problems can arise, however, when configurations with large forces due to short distance between atoms are included in the data set. Here, a significant improvement to the Gaussian process regression approach is obtained by basing the difference measure between two atomic configurations in the covariance function on the inverted inter-atomic distances and by adding a new early stopping criterion for the path relaxation phase. This greatly improves the performance of the method in two applications where the original formulation does not work well: a dissociative adsorption of an H 2 molecule on a Cu(110) Surface and a diffusion hop of an H 2 O molecule on an ice Ih(0001) surface. Also, the revised method works better in the previously analyzed benchmark application to rearrangement transitions of a heptamer island on a surface, requiring fewer energy and force evaluations for convergence to the minimum energy path.

Methodology for finding optimal tunneling paths and evaluating tunneling rates for atomic rearrangements is described. First, an optimal JWKB tunneling path for a system with fixed energy is obtained using a line integral extension of the nudged elastic band method. Then, a calculation of the dynamics along the path is used to determine the temperature at which it corresponds to an optimal Feynman path for thermally activated tunneling (instanton) and a harmonic approximation is used to estimate the transition rate. The method is illustrated with calculations for a modified two-dimensional Müller-Brown surface but is efficient enough to be used in combination with electronic structure calculations of the energy and atomic forces in systems containing many atoms. An example is presented where tunneling is the dominant mechanism well above room temperature as an H3BNH3 molecule dissociates to form H2. Also, a solid-state example is presented where density functional theory calculations of H atom tunneling in a Ta crystal give close agreement with experimental measurements on hydrogen diffusion over a wide range in temperature.

Minimum energy paths for transitions such as atomic and/or spin rearrangements in thermalized systems are the transition paths of largest statistical weight. Such paths are frequently calculated using the nudged elastic band method, where an initial path is iteratively shifted to the nearest minimum energy path. The computational effort can be large, especially when ab initio or electron density functional calculations are used to evaluate the energy and atomic forces. Here, we show how the number of such evaluations can be reduced by an order of magnitude using a Gaussian process regression approach where an approximate energy surface is generated and refined in each iteration. When the goal is to evaluate the transition rate within harmonic transition state theory, the evaluation of the Hessian matrix at the initial and final state minima can be carried out beforehand and used as input in the minimum energy path calculation, thereby improving stability and reducing the number of iterations needed for convergence. A Gaussian process model also provides an uncertainty estimate for the approximate energy surface, and this can be used to focus the calculations on the lesser-known part of the path, thereby reducing the number of needed energy and force evaluations to a half in the present calculations. The methodology is illustrated using the two-dimensional Müller-Brown potential surface and performance assessed on an established benchmark involving 13 rearrangement transitions of a heptamer island on a solid surface.

We introduce a fully stand-alone version of the Quantum Chemistry Electron Ionization Mass Spectra (QCEIMS) program [S. Grimme, Angew. Chem. Int. Ed., 2013, 52, 6306] allowing efficient simulations for molecules composed of elements with atomic numbers up to Z = 86. The recently developed extended tight-binding semi-empirical method GFN-xTB has been combined with QCEIMS, thereby eliminating dependencies on third-party electronic structure software. Furthermore, for reasonable calculations of ionization potentials, as required by the method, a second tight-binding variant, IPEA-xTB, is introduced here. This novel combination of methods allows the automatic, fast and reasonably accurate computation of electron ionization mass spectra for structurally different molecules across the periodic table. In order to validate and inspect the transferability of the method, we perform large-scale simulations for some representative organic, organometallic, and main-group inorganic systems. Theoretical spectra for 23 molecules are compared directly to experimental data taken from standard databases. For the first time, realistic quantum chemistry based EI-MS for organometallic systems like ferrocene or copper(II) acetylacetonate are presented. Compared to previously used semiempirical methods, GFN-xTB is faster, more robust, and yields overall higher quality spectra. The partially analysed theoretical reaction and fragmentation mechanisms are chemically reasonable and reveal in unprecedented detail the extreme complexity of high energy gas phase ion chemistry including complicated rearrangement reactions prior to dissociation.

Hydrogen (H) atom diffusion on dust grain surfaces is the rate limiting step in many hydrogenation reactions taking place in interstellar clouds. In cold (10-30 K) molecular clouds, the dust grains are coated by amorphous water ice. Therefore, H adatom mobility on ice surfaces is of fundamental importance in this context. We have calculated H atom adsorption and diffusion on both crystalline and amorphous ice surfaces using an analytic interaction potential for H2O H. Tunneling rates for H atom hops between adsorption sites are explicitly calculated, the kinetic Monte Carlo method is used to simulate long time scale evolution and the diffusion coefficient, D, is evaluated for the temperature range 5-120 K. For ice Ih, we find D = 1.6 X 10(-7) cm(2)/s at 10 K and below that temperature tunneling becomes the dominant diffusion mechanism. On the amorphous ice surface, the mobility of H is much slower than for ice Ih, D = 5.8 X 10(-11) cm(2)/s at 25 K. Below 25 K, the H adatom becomes trapped in the deepest adsorption site which is located in a small surface pore. Furthermore, by blocking this site the diffusion increases by several orders of magnitude. H, formation is thus likely to take place in deep adsorption sites at the low coverage and low temperature characteristic of molecular clouds. Deep adsorption sites can also explain experimental observations indicating that tunneling does not significantly contribute to the diffusivity, since even though H adatoms can tunnel between shallow adsorption sites, the rate-determining transitions out of deep sites require thermal activation. Furthermore, in the context of coarse-grained astrochemical models, we find the ratio of the activation energy of diffusion and the adsorption energy of H to be 0.64 on amorphous water ice.

We present negative ion-mode simulations within the QCEIMS program [Grimme, Angew. Chem., Int. Ed., 2013, 52, 6306]. It is an exhaustive and robust ab initio molecular dynamics/stochastic algorithm used to perform simulations of unimolecular decomposition of anions, in unprecedented detail. The objective of this approach is to compliment electron attachment spectroscopy and aid in the interpretation of relevant dissociation dynamics. Prototypical simulations are performed for the four nitrile compounds acetonitrile, cyanamide, aminoacetonitrile, and trifluoroacetonitrile. The unique decomposition pathways which naturally occur in the simulations are addressed along with fractional yields, reaction times and relative intensities of the fragments. Furthermore, trajectories of selected decomposition pathways of the aminoacetonitrile anion are investigated in greater detail, where we find that the relevant HOMO of the anion has a mixed π* and σ* character delocalized over the entire molecule.