Exploring the frontiers of condensed-phase chemistry with a general reactive machine learning potential.

Atomistic simulation has a broad range of applications from drug design to materials discovery. Machine learning interatomic potentials (MLIPs) have become an efficient alternative to computationally expensive ab initio simulations. For this reason, chemistry and materials science would greatly benefit from a general reactive MLIP, that is, an MLIP that is applicable to a broad range of reactive chemistry without the need for refitting. Here we develop a general reactive MLIP (ANI-1xnr) through automated sampling of condensed-phase reactions. ANI-1xnr is then applied to study five distinct systems: carbon solid-phase nucleation, graphene ring formation from acetylene, biofuel additives, combustion of methane and the spontaneous formation of glycine from early earth small molecules. In all studies, ANI-1xnr closely matches experiment (when available) and/or previous studies using traditional model chemistry methods. As such, ANI-1xnr proves to be a highly general reactive MLIP for C, H, N and O elements in the condensed phase, enabling high-throughput in silico reactive chemistry experimentation.

Bayesian-Inference-Driven Model Parametrization and Model Selection for 2CLJQ Fluid Models.

A high level of physical detail in a molecular model improves its ability to perform high accuracy simulations but can also significantly affect its complexity and computational cost. In some situations, it is worthwhile to add complexity to a model to capture properties of interest; in others, additional complexity is unnecessary and can make simulations computationally infeasible. In this work, we demonstrate the use of Bayesian inference for molecular model selection, using Monte Carlo sampling techniques accelerated with surrogate modeling to evaluate the Bayes factor evidence for different levels of complexity in the two-centered Lennard-Jones + quadrupole (2CLJQ) fluid model. Examining three nested levels of model complexity, we demonstrate that the use of variable quadrupole and bond length parameters in this model framework is justified only for some chemistries. Through this process, we also get detailed information about the distributions and correlation of parameter values, enabling improved parametrization and parameter analysis. We also show how the choice of parameter priors, which encode previous model knowledge, can have substantial effects on the selection of models, penalizing careless introduction of additional complexity. We detail the computational techniques used in this analysis, providing a roadmap for future applications of molecular model selection via Bayesian inference and surrogate modeling.

Predicting phosphorescence energies and inferring wavefunction localization with machine learning.

Phosphorescence is commonly utilized for applications including light-emitting diodes and photovoltaics. Machine learning (ML) approaches trained on ab initio datasets of singlet-triplet energy gaps may expedite the discovery of phosphorescent compounds with the desired emission energies. However, we show that standard ML approaches for modeling potential energy surfaces inaccurately predict singlet-triplet energy gaps due to the failure to account for spatial localities of spin transitions. To solve this, we introduce localization layers in a neural network model that weight atomic contributions to the energy gap, thereby allowing the model to isolate the most determinative chemical environments. Trained on the singlet-triplet energy gaps of organic molecules, we apply our method to an out-of-sample test set of large phosphorescent compounds and demonstrate the substantial improvement that localization layers have on predicting their phosphorescence energies. Remarkably, the inferred localization weights have a strong relationship with the ab initio spin density of the singlet-triplet transition, and thus infer localities of the molecule that determine the spin transition, despite the fact that no direct electronic information was provided during training. The use of localization layers is expected to improve the modeling of many localized, non-extensive phenomena and could be implemented in any atom-centered neural network model.

Molecular Calculation of the Critical Parameters of Classical Helium.

We compute the vapor-liquid critical coordinates of a model of helium in which nuclear quantum effects are absent. We employ highly accurate ab initio pair and three-body potentials and calculate the critical parameters rigorously in two ways. First, we calculate the virial coefficients up to the seventh and find the point where an isotherm satisfies the critical conditions. Second, we use Gibbs Ensemble Monte Carlo (GEMC) to calculate the vapor-liquid equilibrium, and extrapolate the phase envelope to the critical point. Both methods yield results that are consistent within their uncertainties. The critical temperature of "classical helium" is 13.0 K (compared to 5.2 K for real helium), the critical pressure is 0.93 MPa, and the critical density is 28.4 mol·L-1, with expanded uncertainties (corresponding to a 95% confidence interval) on the order of 0.1 K, 0.02 MPa, and 0.5 mol·L-1, respectively. The effect of three-body interactions on the location of the critical point is small (lowering the critical temperature by roughly 0.1 K), suggesting that we are justified in ignoring four-body and higher interactions in our calculations. This work is motivated by the use of corresponding-states models for mixtures containing helium (such as some natural gases) at higher temperatures where quantum effects are expected to be negligible; in these situations, the distortion of the critical properties by quantum effects causes problems for the corresponding-states treatment.

Uncertainty quantification confirms unreliable extrapolation toward high pressures for united-atom Mie λ-6 force field.

Molecular simulation results at extreme temperatures and pressures can supplement experimental data when developing fundamental equations of state. Since most force fields are optimized to agree with vapor-liquid equilibria (VLE) properties, however, the reliability of the molecular simulation results depends on the validity/transferability of the force field at higher temperatures and pressures. As demonstrated in this study, although state-of-the-art united-atom Mie λ-6 potentials for normal and branched alkanes provide accurate estimates for VLE, they tend to over-predict pressures for dense supercritical fluids and compressed liquids. The physical explanation for this observation is that the repulsive barrier is too steep for the "optimal" united-atom Mie λ-6 potential parameterized with VLE properties. Bayesian inference confirms that no feasible combination of non-bonded parameters (ϵ, σ, and λ) is capable of simultaneously predicting saturated vapor pressures, saturated liquid densities, and pressures at high temperatures and densities. This conclusion has both practical and theoretical ramifications, as more realistic non-bonded potentials may be required for accurate extrapolation to high pressures of industrial interest.

Configuration-Sampling-Based Surrogate Models for Rapid Parameterization of Non-Bonded Interactions.

In this study, we present an approach for rapid force field parameterization and uncertainty quantification of the non-bonded interaction parameters for classical force fields. The accuracy of most thermophysical properties, and especially vapor-liquid equilibria (VLE), obtained from molecular simulation depends strongly on the non-bonded interactions. Traditionally, non-bonded interactions are parameterized to agree with macroscopic properties by performing large amounts of direct molecular simulation. Due to the computational cost of molecular simulation, surrogate models (i.e., efficient models that approximate direct molecular simulation results) are an essential tool for high-dimensional parameterization and uncertainty quantification of non-bonded interactions. The present study compares two different configuration-sampling-based surrogate models, namely, Multistate Bennett Acceptance Ratio (MBAR) and Pair Correlation Function Rescaling (PCFR). MBAR and PCFR are coupled with the Isothermal Isochoric (ITIC) thermodynamic integration method for estimating vapor-liquid saturation properties. We find that MBAR and PCFR are complementary in their roles. Specifically, PCFR is preferred when exploring distant regions of the parameter space while MBAR is better in the local domain.

Uncertainty quantification and propagation of errors of the Lennard-Jones 12-6 parameters for n-alkanes.

Molecular simulation has the ability to predict various physical properties that are difficult to obtain experimentally. For example, we implement molecular simulation to predict the critical constants (i.e., critical temperature, critical density, critical pressure, and critical compressibility factor) for large n-alkanes that thermally decompose experimentally (as large as C48). Historically, molecular simulation has been viewed as a tool that is limited to providing qualitative insight. One key reason for this perceived weakness in molecular simulation is the difficulty to quantify the uncertainty in the results. This is because molecular simulations have many sources of uncertainty that propagate and are difficult to quantify. We investigate one of the most important sources of uncertainty, namely, the intermolecular force field parameters. Specifically, we quantify the uncertainty in the Lennard-Jones (LJ) 12-6 parameters for the CH4, CH3, and CH2 united-atom interaction sites. We then demonstrate how the uncertainties in the parameters lead to uncertainties in the saturated liquid density and critical constant values obtained from Gibbs Ensemble Monte Carlo simulation. Our results suggest that the uncertainties attributed to the LJ 12-6 parameters are small enough that quantitatively useful estimates of the saturated liquid density and the critical constants can be obtained from molecular simulation.

An improved statistical analysis for predicting the critical temperature and critical density with Gibbs ensemble Monte Carlo simulation.

A rigorous statistical analysis is presented for Gibbs ensemble Monte Carlo simulations. This analysis reduces the uncertainty in the critical point estimate when compared with traditional methods found in the literature. Two different improvements are recommended due to the following results. First, the traditional propagation of error approach for estimating the standard deviations used in regression improperly weighs the terms in the objective function due to the inherent interdependence of the vapor and liquid densities. For this reason, an error model is developed to predict the standard deviations. Second, and most importantly, a rigorous algorithm for nonlinear regression is compared to the traditional approach of linearizing the equations and propagating the error in the slope and the intercept. The traditional regression approach can yield nonphysical confidence intervals for the critical constants. By contrast, the rigorous algorithm restricts the confidence regions to values that are physically sensible. To demonstrate the effect of these conclusions, a case study is performed to enhance the reliability of molecular simulations to resolve the n-alkane family trend for the critical temperature and critical density.