DOI: 10.6084/m9.figshare.29949230.v1
title: Landscape17

description: OverviewMachine learning interatomic potentials (MLIPs) have achieved remarkable accuracy on standard benchmarks, yet their ability to reproduce molecular kinetics – critical for reaction rate calculations – remains largely unexplored. We introduce the Landscape17 dataset, which systematically expands upon the rMD17 dataset by providing complete kinetic transition networks (KTNs) for the six molecules within rMD17 that have more than one distinct local minimum structure. This dataset features global potential energy surface representations generated using the energy landscape framework and includes regions crucial for accurately reproducing both thermodynamic and kinetic properties. For each of the selected six molecules (ethanol, malonaldehyde, paracetamol, salicylic acid, azobenzene, and aspirin) we provide all the minima and transition states, along with configurations from the two approximate steepest-descent paths connecting each transition state to the corresponding minima, computed using hybrid-level density functional theory. These paths underpin the most probable routes between minima at finite temperature, offering essential configurations for understanding system kinetics.We utilized TopSearch, an open-source Python package, to perform landscape exploration, at an estimated cost of 10<sup>5 </sup>CPUh. For azobenzene we employed the OPTIM package Cambridge Energy Landscapes software suite. For each KTN, we excluded repeated permutational isomers and structures related by the inversion operation, as these can be reconstructed through symmetry operations.<br>*Note: Hessian eigenvalues should be multiplied by 9442.713128494901 to obtain values in eV/A<sup>2</sup>.PublicationFurther details of this dataset, and associated references, are given in the corresponding manuscript:Cărare, V., Thiemann, F.L., Morrow, J.D. <i>et al.</i> Global properties of the energy landscape: a testing and training arena for machine learned potentials. <i>npj Comput Mater</i> (2025). https://doi.org/10.1038/s41524-025-01878-xIf using this dataset, please consider citing the above.MethodsEnergy landscape frameworkThe energy landscape framework provides a comprehensive approach to mapping surface topography through the identification and characterization of stationary points. These are atomic configurations at which the gradient vanishes and we focus on local minima and the transition states that connect them, which are distinguished by their Hessian eigenvalue spectrum. Local minima exhibit only positive and zero Hessian eigenvalues, indicating that any displacement of internal coordinates increases the energy. Transition states are defined as first-order saddle points with exactly one negative eigenvalue, corresponding to a local maximum along the reaction coordinate, with positive curvature in the orthogonal eigendirections (aside from the zero eigenvalues corresponding to overall rotation and translation). These stationary points can be represented by weighted graphs, known as kinetic transition networks (KTNs), where minima serve as nodes and edges connect minima that are directly linked by transition states. Appropriate post-processing using standard tools of statistical mechanics and unimolecular rate theory enables efficient computation of observable thermodynamic and kinetic properties within well-defined approximations. In particular, the explicit inclusion of transition states, which are more difficult to characterise using standard molecular dynamics, allows for assessment of global kinetics and comparison of MLIP landscapes with the DFT reference.Density functional theory calculationsThe reference potential energy landscapes were computed using density functional theory with the ωB97x hybrid-energy exchange correlation functional and a 6-31G(d) basis set within Psi4. These settings are consistent with the ones used to generate the ANI2x training data. We applied tight energy and density convergence criteria (E CONVERGENCE and D CONVERGENCE) of 10<sup>−9</sup> Hartree and 10<sup>−9</sup> a.u., along with extremely fine integration grids (100 radial and 770 spherical points) and a restricted Kohn-Sham reference. We validated the convergence, grid, and spin settings against published data from rMD17, using the appropriate functional and basis set: PBE/def2-SVP. We achieved energies and forces within 0.1 meV/atom and 5 meV/Å respectively, well within the standard acceptable resolution of 1 meV/atom and 10 meV/Å.Minima identificationLocal minima were first collected from the basin-hopping global optimization runs, an approach that has has been employed successfully for diverse molecular and abstract landscapes. The basin-hopping exploration employed random angular perturbations applied to flexible dihedral angles, as identified by the Atomic Simulation Environment package. These surveys used 100 basin-hopping steps with an accept/reject Metropolis condition equivalent to a temperature of 100 K. The convergence criteria required either the maximum force component to fall below 10<sup>−3</sup> eV/Å or the relative energy change between steps to drop below 10<sup>7</sup> eV multiplied by machine precision.Transition state locationTransition state searches were performed between each minimum and its three nearest neighbours, determined from the Euclidean distance with optimal alignment via the MINPERMDIST routine. This procedure minimises the distance with respect to translation, rotation, permutation and, additionally, we include the inversion operation. The alignment is not deterministic when permutations are included, and employs a shortest augmenting path algorithm inside an iterative loop. Distinct stationary points were distinguished by a root mean square distance threshold greater than 0.3 Å and energy differences exceeding 10<sup>−3</sup> eV. Only one representative permutation-inversion isomer was retained for each minimum and transition state. Transition states were located using a two-step protocol starting with a nudged elastic band (NEB) calculation. Initial pathways between minima were generated through linear interpolation in internal coordinates after endpoint alignment. The NEB algorithm optimizes these interpolations of 20 images with a 50 eV/Å spring constant and convergence criterion of 10<sup>−2</sup> eV/Å for the maximum force component. This double-ended phase of the calculation only needs to converge sufficiently to identify the local maxima in the profile, which are taken as starting points for accurate refinement using hybrid eigenvector-following. Here the smallest non-zero eigenvalue and the corresponding eigenvector are obtained using a variational method, with minimisation in all orthogonal directions. Convergence required RMS forces below 3×10<sup>−2</sup> eV/Å. Transition states were verified through Hessian analysis, confirming exactly one negative eigenvalue. Only transition states with barriers below 1 eV were retained. We do not assume that a path links the two original end minima from the NEB phase, even when there is only a single transition state in the profile. The connectivity of each transition state is always established from the corresponding pathways. In general, there may be multiple transition states between two minima, and there may be gaps in the connection profile. The two-phase procedure is applied until a complete discrete path is obtained, using the missing connection algorithm to propose new pairs of minima for additional searches.Pathways and network constructionEach transition state connects two local minima, identified through minimization from perturbed transition state geometries. Perturbations of approximately 0.3<b>y</b> (0.6<b>y</b> for flatter modes such as the cis-cis transition state in azobenzene) were applied along parallel and antiparallel to the normalised eigenvector, <b>y</b>, corresponding to the negative eigenvalue, followed by LBFGS minimization. All intermediate configurations during minimization were stored to map the approximate steepest-descent paths from transition states to connected minima, excluding initial configurations with force components exceeding 2 eV/Å.