A fragmented, permutationally invariant polynomial approach for potential energy surfaces of large molecules: Application to N-methyl acetamide

Extending atomic decomposition and many-body representation with a chemistry-motivated approach to machine learning potentials

Most widely used machine learning potentials for condensed-phase applications rely on many-body permutationally invariant polynomial or atom-centered neural networks. However, these approaches face challenges in achieving chemical interpretability in atomistic energy decomposition and fully matching the computational efficiency of traditional force fields. Here we present a method that combines aspects of both approaches and balances accuracy and force-field-level speed. This method utilizes a monomer-centered representation, where the potential energy is decomposed into the sum of chemically meaningful monomeric energies. The structural descriptors of monomers are described by one-body and two-body effective interactions, enforced by appropriate sets of permutationally invariant polynomials as inputs to the feed-forward neural networks. Systematic assessments of models for gas-phase water trimer, liquid water, methane-water cluster and liquid carbon dioxide are performed. The improved accuracy, efficiency and flexibility of this method have promise for constructing accurate machine learning potentials and enabling large-scale quantum and classical simulations for complex molecular systems.

MOLPIPx: An end-to-end differentiable package for permutationally invariant polynomials in Python and Rust

In this work, we present MOLPIPx, a versatile library designed to seamlessly integrate permutationally invariant polynomials with modern machine learning frameworks, enabling the efficient development of linear models, neural networks, and Gaussian process models. These methodologies are widely employed for parameterizing potential energy surfaces across diverse molecular systems. MOLPIPx leverages two powerful automatic differentiation engines-JAX and EnzymeAD-Rust-to facilitate the efficient computation of energy gradients and higher-order derivatives, which are essential for tasks such as force field development and dynamic simulations. MOLPIPx is available at https://github.com/ChemAI-Lab/molpipx.

New Algorithms to Generate Permutationally Invariant Polynomials and Fundamental Invariants for Potential Energy Surface Fitting

Symmetric functions, such as Permutationally Invariant Polynomials (PIPs) and Fundamental Invariants (FIs), are effective and concise descriptors for incorporating permutation symmetry into neural network (NN) potential energy surface (PES) fitting. The traditional algorithm for generating such symmetric polynomials has a factorial time complexity of N!, where N is the number of identical atoms, posing a significant challenge to applying symmetric polynomials as descriptors of NN PESs for larger systems, particularly with more than 10 atoms. Herein, we report a new algorithm which has only linear time complexity for identical atoms. It can tremendously accelerate generation process of symmetric polynomials for molecular systems. The proposed algorithm is based on graph connectivity analysis following the action of the generation set of molecular permutational group. For instance, in the case of calculating the invariant polynomials for a 15-atom molecule, such as tropolone, our algorithm is approximately 2 million times faster than the previous method. The efficiency of the new algorithm can be further enhanced with increasing molecular size and number of identical atoms, making the FI-NN approach feasible for systems with over 10 atoms and high symmetry demands.

DFT-Based Permutationally Invariant Polynomial Potentials Capture the Twists and Turns of C14H30

Hydrocarbons are ubiquitous as fuels, solvents, lubricants, and as the principal components of plastics and fibers, yet our ability to predict their dynamical properties is limited to force-field mechanics. Here, we report two machine-learned potential energy surfaces (PESs) for the linear 44-atom hydrocarbon C14H30 using an extensive data set of roughly 250,000 density functional theory (DFT) (B3LYP) energies for a large variety of configurations, obtained using MM3 direct-dynamics calculations at 500, 1000, and 2500 K. The surfaces, based on Permutationally Invariant Polynomials (PIPs) and using both a many-body expansion approach and a fragmented-basis approach, produce precise fits for energies and forces and also produce excellent out-of-sample agreement with direct DFT calculations for torsional and dihedral angle potentials. Going beyond precision, the PESs are used in molecular dynamics calculations that demonstrate the robustness of the PESs for a large range of conformations. The many-body PIPs PES, although more compute intensive than the fragmented-basis one, is directly transferable for other linear hydrocarbons.

Machine Learning of Reactive Potentials

In the past two decades, machine learning potentials (MLPs) have driven significant developments in chemical, biological, and material sciences. The construction and training of MLPs enable fast and accurate simulations and analysis of thermodynamic and kinetic properties. This review focuses on the application of MLPs to reaction systems with consideration of bond breaking and formation. We review the development of MLP models, primarily with neural network and kernel-based algorithms, and recent applications of reactive MLPs (RMLPs) to systems at different scales. We show how RMLPs are constructed, how they speed up the calculation of reactive dynamics, and how they facilitate the study of reaction trajectories, reaction rates, free energy calculations, and many other calculations. Different data sampling strategies applied in building RMLPs are also discussed with a focus on how to collect structures for rare events and how to further improve their performance with active learning.

No Headache for PIPs: A PIP Potential for Aspirin Runs Much Faster and with Similar Precision Than Other Machine-Learned Potentials

Assessments of machine-learning (ML) potentials are an important aspect of the rapid development of this field. We recently reported an assessment of the linear-regression permutationally invariant polynomial (PIP) method for ethanol, using the widely used (revised) rMD17 data set. We demonstrated that the PIP approach outperformed numerous other methods, e.g., ANI, PhysNet, sGDML, and p-KRR, with respect to precision and notably with respect to speed [Houston et al., J. Chem. Phys. 2022, 156, 044120]. Here, we extend this assessment to the 21-atom aspirin molecule, using the rMD17 data set, with a focus on the speed of evaluation. Both energies and forces are used for training, and the precision of several PIPs is examined for both. Normal mode frequencies, the methyl torsional potential, and 1d vibrational energies for an OH stretch are presented. We show that the PIP approach achieves the level of precision obtained from other ML methods, e.g., atom-centered neural network methods, linear regression ACE, and kernel methods, as reported by Kovács et al. in J. Chem. Theory Comput. 2021, 17, 7696-7711. More significantly, we show that the PIP PESs run much faster than all other ML methods, whose timings were evaluated in that paper. We also show that the PIP PES extrapolates well enough to describe several internal motions of aspirin, including an OH stretch.

Accurate fundamental invariant-neural network representation of ab initio potential energy surfaces

Highly accurate potential energy surfaces are critically important for chemical reaction dynamics. The large number of degrees of freedom and the intricate symmetry adaption pose a big challenge to accurately representing potential energy surfaces (PESs) for polyatomic reactions. Recently, our group has made substantial progress in this direction by developing the fundamental invariant-neural network (FI-NN) approach. Here, we review these advances, demonstrating that the FI-NN approach can represent highly accurate, global, full-dimensional PESs for reactive systems with even more than 10 atoms. These multi-channel reactions typically involve many intermediates, transition states, and products. The complexity and ruggedness of this potential energy landscape present even greater challenges for full-dimensional PES representation. These PESs exhibit a high level of complexity, molecular size, and accuracy of fit. Dynamics simulations based on these PESs have unveiled intriguing and novel reaction mechanisms, providing deep insights into the intricate dynamics involved in combustion, atmospheric, and organic chemistry.

MBX: A many-body energy and force calculator for data-driven many-body simulations

Many-Body eXpansion (MBX) is a C++ library that implements many-body potential energy functions (PEFs) within the "many-body energy" (MB-nrg) formalism. MB-nrg PEFs integrate an underlying polarizable model with explicit machine-learned representations of many-body interactions to achieve chemical accuracy from the gas to the condensed phases. MBX can be employed either as a stand-alone package or as an energy/force engine that can be integrated with generic software for molecular dynamics and Monte Carlo simulations. MBX is parallelized internally using Open Multi-Processing and can utilize Message Passing Interface when available in interfaced molecular simulation software. MBX enables classical and quantum molecular simulations with MB-nrg PEFs, as well as hybrid simulations that combine conventional force fields and MB-nrg PEFs, for diverse systems ranging from small gas-phase clusters to aqueous solutions and molecular fluids to biomolecular systems and metal-organic frameworks.

Accurate and Interpretable Dipole Interaction Model-Based Machine Learning for Molecular Polarizability

Polarizabilities play significant roles in describing dispersive and inductive interactions of the atom and molecular systems. However, an accurate prediction of molecular polarizabilities from first principles is computationally prohibitive. Although physical models or statistical machine learning models have been proposed, either a lack of accurate description of local chemical environments or demanding a large number of samples for training has limited their practical applications. In this study, we combine a physically inspired dipole interaction model and an accurate neural network method for predicting the polarizability tensors of molecules. With the local chemical environment precisely described and the requirement of rotational covariance naturally fulfilled, this hybrid model is proven to give an accurate molecular polarizability prediction, essentially reducing the number of training samples. The atomic polarizabilities are physically interpretable and transferable to larger molecules unseen in the training set. This promising method may find its wide range of applications, such as spectroscopic simulations and the construction of polarizable force fields.

PESPIP: Software to fit complex molecular and many-body potential energy surfaces with permutationally invariant polynomials

We wish to describe a potential energy surface by using a basis of permutationally invariant polynomials whose coefficients will be determined by numerical regression so as to smoothly fit a dataset of electronic energies as well as, perhaps, gradients. The polynomials will be powers of transformed internuclear distances, usually either Morse variables, exp(-r_i,j/λ), where λ is a constant range hyperparameter, or reciprocals of the distances, 1/r_i,j. The question we address is how to create the most efficient basis, including (a) which polynomials to keep or discard, (b) how many polynomials will be needed, (c) how to make sure the polynomials correctly reproduce the zero interaction at a large distance, (d) how to ensure special symmetries, and (e) how to calculate gradients efficiently. This article discusses how these questions can be answered by using a set of programs to choose and manipulate the polynomials as well as to write efficient Fortran programs for the calculation of energies and gradients. A user-friendly interface for access to monomial symmetrization approach results is also described. The software for these programs is now publicly available.

Δ-Machine Learned Potential Energy Surfaces and Force Fields

There has been great progress in developing machine-learned potential energy surfaces (PESs) for molecules and clusters with more than 10 atoms. Unfortunately, this number of atoms generally limits the level of electronic structure theory to less than the "gold standard" CCSD(T) level. Indeed, for the well-known MD17 dataset for molecules with 9-20 atoms, all of the energies and forces were obtained with DFT calculations (PBE). This Perspective is focused on a Δ-machine learning method that we recently proposed and applied to bring DFT-based PESs to close to CCSD(T) accuracy. This is demonstrated for hydronium, N-methylacetamide, acetyl acetone, and ethanol. For 15-atom tropolone, it appears that special approaches (e.g., molecular tailoring, local CCSD(T)) are needed to obtain the CCSD(T) energies. A new aspect of this approach is the extension of Δ-machine learning to force fields. The approach is based on many-body corrections to polarizable force field potentials. This is examined in detail using the TTM2.1 water potential. The corrections make use of our recent CCSD(T) datasets for 2-b, 3-b, and 4-b interactions for water. These datasets were used to develop a new fully ab initio potential for water, termed q-AQUA.

Hierarchical Materials from High Information Content Macromolecular Building Blocks: Construction, Dynamic Interventions, and Prediction

Hierarchical materials that exhibit order over multiple length scales are ubiquitous in nature. Because hierarchy gives rise to unique properties and functions, many have sought inspiration from nature when designing and fabricating hierarchical matter. More and more, however, nature's own high-information content building blocks, proteins, peptides, and peptidomimetics, are being coopted to build hierarchy because the information that determines structure, function, and interfacial interactions can be readily encoded in these versatile macromolecules. Here, we take stock of recent progress in the rational design and characterization of hierarchical materials produced from high-information content blocks with a focus on stimuli-responsive and "smart" architectures. We also review advances in the use of computational simulations and data-driven predictions to shed light on how the side chain chemistry and conformational flexibility of macromolecular blocks drive the emergence of order and the acquisition of hierarchy and also on how ionic, solvent, and surface effects influence the outcomes of assembly. Continued progress in the above areas will ultimately usher in an era where an understanding of designed interactions, surface effects, and solution conditions can be harnessed to achieve predictive materials synthesis across scale and drive emergent phenomena in the self-assembly and reconfiguration of high-information content building blocks.

Semiclassical and VSCF/VCI Calculations of the Vibrational Energies of trans- and gauche-Ethanol Using a CCSD(T) Potential Energy Surface

A recent full-dimensional Δ-Machine learning potential energy surface (PES) for ethanol is employed in semiclassical and vibrational self-consistent field (VSCF) and virtual-state configuration interaction (VCI) calculations, using MULTIMODE, to determine the anharmonic vibrational frequencies of vibration for both the trans and gauche conformers of ethanol. Both semiclassical and VSCF/VCI energies agree well with the experimental data. We find significant mixing between the VSCF basis states due to Fermi resonances between bending and stretching modes. The same effects are also accurately described by the full-dimensional semiclassical calculations. These are the first high-level anharmonic calculations using a PES, in particular a "gold-standard" CCSD(T) one.

Data-Driven Many-Body Potential Energy Functions for Generic Molecules: Linear Alkanes as a Proof-of-Concept Application

We present a generalization of the many-body energy (MB-nrg) theoretical/computational framework that enables the development of data-driven potential energy functions (PEFs) for generic covalently bonded molecules, with arbitrary quantum mechanical accuracy. The "nearsightedness of electronic matter" is exploited to define monomers as "natural building blocks" on the basis of their distinct chemical identity. The energy of generic molecules is then expressed as a sum of individual many-body energies of incrementally larger subsystems. The MB-nrg PEFs represent the low-order n-body energies, with n = 1-4, using permutationally invariant polynomials derived from electronic structure data carried out at an arbitrary quantum mechanical level of theory, while all higher-order n-body terms (n > 4) are represented by a classical many-body polarization term. As a proof-of-concept application of the general MB-nrg framework, we present MB-nrg PEFs for linear alkanes. The MB-nrg PEFs are shown to accurately reproduce reference energies, harmonic frequencies, and potential energy scans of alkanes, independently of their length. Since, by construction, the MB-nrg framework introduced here can be applied to generic covalently bonded molecules, we envision future computer simulations of complex molecular systems using data-driven MB-nrg PEFs, with arbitrary quantum mechanical accuracy.

A CCSD(T)-Based 4-Body Potential for Water

High-level, ab initio calculations find that the 4-body (4-b) interaction is needed to account for near-100% of the total interaction energy for water clusters as large as the 21-mer. Motivated by this, we report a permutationally invariant polynomial potential energy surface (PES) for the 4-body interaction. This machine-learned PES is a fit to 2119 symmetry-unique, CCSD(T)-F12a/haTZ 4-b interaction energies. Configurations for these come from tetramer direct-dynamics calculations, fragments from an MD water simulation at 300 K, and tetramer fragments in a variety of water clusters. The PIP basis is purified to ensure that the PES goes rigorously to zero in monomer+trimer and dimer+dimer dissociations. The 4-b energies of isomers of the hexamer calculated with the new PES are shown to be in better agreement with benchmark CCSD(T) results than those from the MB-pol potential. Tests on larger clusters further validate the high-fidelity of the PES. The PES is shown to be fast to evaluate, taking 2.4 s for 10⁵ evaluations on a single core of 2.4 GHz Intel Xeon processor, and significantly faster using a parallel version of the PES.

Parsimonious Potential Energy Surface Expansions Using Dictionary Learning with Multipass Greedy Selection

Potential energy surfaces fit with basis set expansions have been shown to provide accurate representations of electronic energies and have enabled a variety of high-accuracy dynamics, kinetics, and spectroscopy applications. The number of terms in these expansions scales poorly with system size, a drawback that challenges their use for systems with more than ∼10 atoms. A solution is presented here using dictionary learning. Subsets of the full set of conventional basis functions are optimized using a newly developed multipass greedy regression method inspired by forward and backward selection methods from the statistics, signal processing, and machine learning literatures. The optimized representations have accuracies comparable to the full set but are 1 or more orders of magnitude smaller, and notably, the number of terms in the optimized multipass greedy expansions scales approximately linearly with the number of atoms.

Breaking the Coupled Cluster Barrier for Machine-Learned Potentials of Large Molecules: The Case of 15-Atom Acetylacetone

Machine-learned potential energy surfaces (PESs) for molecules with more than 10 atoms are typically forced to use lower-level electronic structure methods such as density functional theory (DFT) and second-order Møller-Plesset perturbation theory (MP2). While these are efficient and realistic, they fall short of the accuracy of the "gold standard" coupled-cluster method, especially with respect to reaction and isomerization barriers. We report a major step forward in applying a Δ-machine learning method to the challenging case of acetylacetone, whose MP2 barrier height for H-atom transfer is low by roughly 1.1 kcal/mol relative to the benchmark CCSD(T) barrier of 3.2 kcal/mol. From a database of 2151 local CCSD(T) energies and training with as few as 430 energies, we obtain a new PES with a barrier of 3.5 kcal/mol in agreement with the LCCSD(T) barrier of 3.5 kcal/mol and close to the benchmark value. Tunneling splittings due to H-atom transfer are calculated using this new PES, providing improved estimates over previous ones obtained using an MP2-based PES.

"Dividing and Conquering" and "Caching" in Molecular Modeling

Molecular modeling is widely utilized in subjects including but not limited to physics, chemistry, biology, materials science and engineering. Impressive progress has been made in development of theories, algorithms and software packages. To divide and conquer, and to cache intermediate results have been long standing principles in development of algorithms. Not surprisingly, most important methodological advancements in more than half century of molecular modeling are various implementations of these two fundamental principles. In the mainstream classical computational molecular science, tremendous efforts have been invested on two lines of algorithm development. The first is coarse graining, which is to represent multiple basic particles in higher resolution modeling as a single larger and softer particle in lower resolution counterpart, with resulting force fields of partial transferability at the expense of some information loss. The second is enhanced sampling, which realizes "dividing and conquering" and/or "caching" in configurational space with focus either on reaction coordinates and collective variables as in metadynamics and related algorithms, or on the transition matrix and state discretization as in Markov state models. For this line of algorithms, spatial resolution is maintained but results are not transferable. Deep learning has been utilized to realize more efficient and accurate ways of "dividing and conquering" and "caching" along these two lines of algorithmic research. We proposed and demonstrated the local free energy landscape approach, a new framework for classical computational molecular science. This framework is based on a third class of algorithm that facilitates molecular modeling through partially transferable in resolution "caching" of distributions for local clusters of molecular degrees of freedom. Differences, connections and potential interactions among these three algorithmic directions are discussed, with the hope to stimulate development of more elegant, efficient and reliable formulations and algorithms for "dividing and conquering" and "caching" in complex molecular systems.

Full-dimensional potential energy surface for acetylacetone and tunneling splittings

New, full-dimensional potential energy surface for acetylacetone allows for description of H-tunneling dynamics and characterization of stationary points.

Advances and New Challenges to Bimolecular Reaction Dynamics Theory

Dynamics of bimolecular reactions in the gas phase are of foundational importance in combustion, atmospheric chemistry, interstellar chemistry, and plasma chemistry. These collision-induced chemical transformations are a sensitive probe of the underlying potential energy surface(s). Despite tremendous progress in past decades, our understanding is still not complete. In this Perspective, we survey the recent advances in theoretical characterization of bimolecular reaction dynamics, stimulated by new experimental observations, and identify key new challenges.

Machine Learning Force Fields and Coarse-Grained Variables in Molecular Dynamics: Application to Materials and Biological Systems

Machine learning encompasses tools and algorithms that are now becoming popular in almost all scientific and technological fields. This is true for molecular dynamics as well, where machine learning offers promises of extracting valuable information from the enormous amounts of data generated by simulation of complex systems. We provide here a review of our current understanding of goals, benefits, and limitations of machine learning techniques for computational studies on atomistic systems, focusing on the construction of empirical force fields from ab initio databases and the determination of reaction coordinates for free energy computation and enhanced sampling.

He Inclusion in Ice-like and Clathrate-like Frameworks: A Benchmark Quantum Chemistry Study of Guest-Host Interactions

Energetics and structural properties of selected type and size He@hydrate frameworks, e.g., from regular structured ice channels to clathrate-like cages, are presented from first-principles quantum chemistry methods. The scarcity of information on He@hydrates makes such complexes challenging targets, while their computational study entails an interesting and arduous task. Some of them have been synthesized in the laboratory, which motivates further investigations on their stability. Hence, the main focus is to examine the performance and accuracy of different wave function-based electronic structure methods, such as MP2, CCSD(T), their explicitly correlated (F12) and domain-based local pair-natural orbital (DLPNO) analogs, as well as modern and conventional density functional theory (DFT) approaches, and analytical model potentials available. Different structures are considered, starting from the "simplest system" formed by a noble gas atom (such as He) and one water molecule, followed by the study of the "fundamental units" present in all ice-like and clathrate-like frameworks (such as pentamers and hexamers) and finally the description of interactions in the "building blocks" of three-dimensional (3D) ice channels (e.g., horizontal and perpendicular ice II and I_h) and clathrate-like cages, such as the 5¹² present in the most common sI, sII, and sH clathrate-hydrate structures. The idea is to provide well-converged DLPNO-CCSD(T) and DFMP2/CBS reference datasets that in turn are used to validate how DFT functionals (in total, 29 approaches from generalized-gradient approximation (GGA), meta-GGA, to hybrid and range-separated functionals, including dispersion correction treatments, were checked) and analytical semiempirical/ab initio-based potentials perform compared with high-level alternatives. Within all tested approaches, those best-performing were identified and classified. Most of the DFT/DFT-D functionals, as well as available analytical pairwise model potentials, face difficulties in describing both hydrogen-bonded water frameworks and dispersion bound He-water interactions. Including dispersion corrections yields an overall well-balanced performance for LCωPBE-D3BJ and PBE0-D4 functionals. Such benchmark datasets can benefit research into the development of new cheminformatics models, as can serve to guide and cross-check methodologies, lending increased predicted power to future molecular simulations for investigating the role of structures and phase transitions from nanoscale clusters to macroscopic crystalline structures.

Fitting potential energy surfaces with fundamental invariant neural network. II. Generating fundamental invariants for molecular systems with up to ten atoms

Symmetry adaptation is crucial in representing a permutationally invariant potential energy surface (PES). Due to the rapid increase in computational time with respect to the molecular size, as well as the reliance on the algebra software, the previous neural network (NN) fitting with inputs of fundamental invariants (FIs) has practical limits. Here, we report an improved and efficient generation scheme of FIs based on the computational invariant theory and parallel program, which can be readily used as the input vector of NNs in fitting high-dimensional PESs with permutation symmetry. The newly developed method significantly reduces the evaluation time of FIs, thereby extending the FI-NN method for constructing highly accurate PESs to larger systems beyond five atoms. Because of the minimum size of invariants used in the inputs of the NN, the NN structure can be very flexible for FI-NN, which leads to small fitting errors. The resulting FI-NN PES is much faster on evaluating than the corresponding permutationally invariant polynomial-NN PES.

Efficient Generation of Permutationally Invariant Potential Energy Surfaces for Large Molecules

An efficient method is described for generating a fragmented, permutationally invariant polynomial basis to fit electronic energies and, if available, gradients for large molecules. The method presented rests on the fragmentation of a large molecule into any number of fragments while maintaining the permutational invariance and uniqueness of the polynomials. The new approach improves on a previous one reported by Qu and Bowman by avoiding repetition of polynomials in the fitting basis set and speeding up gradient evaluations while keeping the accuracy of the PES. The method is demonstrated for CH₃–NH–CO–CH₃ (N-methylacetamide) and NH₂–CH₂–COOH (glycine).

From ab initio data to high-dimensional potential energy surfaces: A critical overview and assessment of the development of permutationally invariant polynomial potential energy surfaces for single molecules

The representation of high-dimensional potential energy surfaces by way of the many-body expansion and permutationally invariant polynomials has become a well-established tool for improving the resolution and extending the scope of molecular simulations. The high level of accuracy that can be attained by these potential energy functions (PEFs) is due in large part to their specificity: for each term in the many-body expansion, a species-specific training set must be generated at the desired level of theory and a number of fits attempted in order to obtain a robust and reliable PEF. In this work, we attempt to characterize the numerical aspects of the fitting problem, addressing questions which are of simultaneous practical and fundamental importance. These include concrete illustrations of the nonconvexity of the problem, the ill-conditionedness of the linear system to be solved and possible need for regularization, the sensitivity of the solutions to the characteristics of the training set, and limitations of the approach with respect to accuracy and the types of molecules that can be treated. In addition, we introduce a general approach to the generation of training set configurations based on the familiar harmonic approximation and evaluate the possible benefits to the use of quasirandom sequences for sampling configuration space in this context. Using sulfate as a case study, the findings are largely generalizable and expected to ultimately facilitate the efficient development of PIP-based many-body PEFs for general systems via automation.