M-Chem: a Modular Software Package for Molecular Simulation that Spans Scientific Domains.

We present a new software package called M-Chem that is designed from scratch in C++ and parallelized on shared-memory multi-core architectures to facilitate efficient molecular simulations. Currently, M-Chem is a fast molecular dynamics (MD) engine that supports the evaluation of energies and forces from two-body to many-body all-atom potentials, reactive force fields, coarse-grained models, combined quantum mechanics molecular mechanics (QM/MM) models, and external force drivers from machine learning, augmented by algorithms that are focused on gains in computational simulation times. M-Chem also includes a range of standard simulation capabilities including thermostats, barostats, multi-timestepping, and periodic cells, as well as newer methods such as fast extended Lagrangians and high quality electrostatic potential generation. At present M-Chem is a developer friendly environment in which we encourage new software contributors from diverse fields to build their algorithms, models, and methods in our modular framework. The long-term objective of M-Chem is to create an interdisciplinary platform for computational methods with applications ranging from biomolecular simulations, reactive chemistry, to materials research.

NewtonNet: a Newtonian message passing network for deep learning of interatomic potentials and forces.

We report a new deep learning message passing network that takes inspiration from Newton's equations of motion to learn interatomic potentials and forces. With the advantage of directional information from trainable force vectors, and physics-infused operators that are inspired by Newtonian physics, the entire model remains rotationally equivariant, and many-body interactions are inferred by more interpretable physical features. We test NewtonNet on the prediction of several reactive and non-reactive high quality ab initio data sets including single small molecules, a large set of chemically diverse molecules, and methane and hydrogen combustion reactions, achieving state-of-the-art test performance on energies and forces with far greater data and computational efficiency than other deep learning models.

A benchmark dataset for Hydrogen Combustion.

The generation of reference data for deep learning models is challenging for reactive systems, and more so for combustion reactions due to the extreme conditions that create radical species and alternative spin states during the combustion process. Here, we extend intrinsic reaction coordinate (IRC) calculations with ab initio MD simulations and normal mode displacement calculations to more extensively cover the potential energy surface for 19 reaction channels for hydrogen combustion. A total of ∼290,000 potential energies and ∼1,270,000 nuclear force vectors are evaluated with a high quality range-separated hybrid density functional, ωB97X-V, to construct the reference data set, including transition state ensembles, for the deep learning models to study hydrogen combustion reaction.

Can electric fields drive chemistry for an aqueous microdroplet?

Reaction rates of common organic reactions have been reported to increase by one to six orders of magnitude in aqueous microdroplets compared to bulk solution, but the reasons for the rate acceleration are poorly understood. Using a coarse-grained electron model that describes structural organization and electron densities for water droplets without the expense of ab initio methods, we investigate the electric field distributions at the air-water interface to understand the origin of surface reactivity. We find that electric field alignments along free O-H bonds at the surface are ~16 MV/cm larger on average than that found for O-H bonds in the interior of the water droplet. Furthermore, electric field distributions can be an order of magnitude larger than the average due to non-linear coupling of intramolecular solvent polarization with intermolecular solvent modes which may contribute to even greater surface reactivity for weakening or breaking chemical bonds at the droplet surface.

Protein C-GeM: A Coarse-Grained Electron Model for Fast and Accurate Protein Electrostatics Prediction.

The electrostatic potential (ESP) is a powerful property for understanding and predicting electrostatic charge distributions that drive interactions between molecules. In this study, we compare various charge partitioning schemes including fitted charges, density-based quantum mechanical (QM) partitioning schemes, charge equilibration methods, and our recently introduced coarse-grained electron model, C-GeM, to describe the ESP for protein systems. When benchmarked against high quality density functional theory calculations of the ESP for tripeptides and the crambin protein, we find that the C-GeM model is of comparable accuracy to ab initio charge partitioning methods, but with orders of magnitude improvement in computational efficiency since it does not require either the electron density or the electrostatic potential as input.

Proton Traffic Jam: Effect of Nanoconfinement and Acid Concentration on Proton Hopping Mechanism.

The properties of the water network in concentrated HCl acid pools in nanometer-sized reverse nonionic micelles were probed with TeraHertz absorption, dielectric relaxation spectroscopy, and reactive force field simulations capable of describing proton hopping mechanisms. We identify that only at a critical micelle size of W0 =9 do solvated proton complexes form in the water pool, accompanied by a change in mechanism from Grotthuss forward shuttling to one that favors local oscillatory hopping. This is due to a preference for H+ and Cl- ions to adsorb to the micelle interface, together with an acid concentration effect that causes a "traffic jam" in which the short-circuiting of the hydrogen-bonding motif of the hydronium ion decreases the forward hopping rate throughout the water interior even as the micelle size increases. These findings have implications for atmospheric chemistry, biochemical and biophysical environments, and energy materials, as transport of protons vital to these processes can be suppressed due to confinement, aggregation, and/or concentration.

Recent Advances for Improving the Accuracy, Transferability, and Efficiency of Reactive Force Fields.

Reactive force fields provide an affordable model for simulating chemical reactions at a fraction of the cost of quantum mechanical approaches. However, classically accounting for chemical reactivity often comes at the expense of accuracy and transferability, while computational cost is still large relative to nonreactive force fields. In this Perspective, we summarize recent efforts for improving the performance of reactive force fields in these three areas with a focus on the ReaxFF theoretical model. To improve accuracy, we describe recent reformulations of charge equilibration schemes to overcome unphysical long-range charge transfer, new ReaxFF models that account for explicit electrons, and corrections for energy conservation issues of the ReaxFF model. To enhance transferability we also highlight new advances to include explicit treatment of electrons in the ReaxFF and hybrid nonreactive/reactive simulations that make it possible to model charge transfer, redox chemistry, and large systems such as reverse micelles within the framework of a reactive force field. To address the computational cost, we review recent work in extended Lagrangian schemes and matrix preconditioners for accelerating the charge equilibration method component of ReaxFF and improvements in its software performance in LAMMPS.

A Reactive Force Field with Coarse-Grained Electrons for Liquid Water.

Nonreactive force fields are defined by perturbations of electron density that are relatively small, whereas chemical reactivity involves wholesale electronic rearrangements that make and break bonds. Thus, reactive force fields are incredibly difficult to develop compared to nonreactive force fields, yet at the same time, they fill a critical need when ab initio molecular dynamics methods are not affordable. We introduce a new reactive force field model for water that combines modified nonbonded terms of the ReaxFF model and its embedding in the electrostatic interactions described by our recently introduced coarse-grained electron model (C-GeM). The ReaxFF/C-GeM force field is characterized for many energetic and dissociative water properties for water clusters, structure, and dynamical properties under ambient conditions in the condensed phase, as well as the temperature dependence of density and water diffusion, with very good agreement with experiment. The ReaxFF/C-GeM force field should be more transferable and more broadly applicable to a range of reactive systems involving both proton and electron transfer in the condensed phase.

Stochastic Constrained Extended System Dynamics for Solving Charge Equilibration Models.

We present a new stochastic extended Lagrangian molecular dynamics solution to charge equilibration that eliminates self-consistent field (SCF) calculations, thus eliminating the computational bottleneck in solving the charge distribution with standard SCF solvers. By formulating both charges and chemical potential as latent variables and introducing a holonomic constraint that satisfies charge conservation, the SC-XLMD method accurately reproduces thermodynamic, dynamic, and structural properties within the framework of ReaxFF for a bulk water system and highly reactive RDX molecules simulated at high temperature. The SC-XLMD method shows excellent computational performance and is available in the publicly available LAMMPS package.

C-GeM: Coarse-Grained Electron Model for Predicting the Electrostatic Potential in Molecules.

We have developed a new coarse-grained electron model, C-GeM, in which atoms are represented by a positive core and an electron shell described by Gaussian charge distributions, with the interaction energy between the core and shell reflecting the electronegativity of a given atomic element. By minimizing the electronic shell positions in the field of atomic core positions, the model can provide accurate electrostatic properties of molecules and their interactions. We have tested the performance of the C-GeM model for a set of molecules containing H, C, O, and Cl atoms to show that it can predict the electrostatic potential with high accuracy, and correctly describe the dissociation of HCl into ionic fragments in solution and to neutral atoms in the gas phase. The resulting C-GeM approach offers many advantages over expensive ab initio methods and reactive force field charge equilibration methodologies: it can rapidly predict the electrostatic potential surfaces of molecules, molecules dissociate into integer charge fragments so that redox reactions are easily described, there is no unphysical long-range charge transfer, it can account for out-of-plane polarization, and charges are not required to be centered on atoms, thereby accounting for electrostatic features such as sigma holes.

Inertial extended-Lagrangian scheme for solving charge equilibration models.

The inertial extended Lagrangian/self-consistent field scheme (iEL-SCF) has been adopted for solving charge equilibration in LAMMPS as part of the reactive force field ReaxFF, which due to the charge conservation constraint requires solving two sets of linear system of equations for the new charges at each molecular dynamics time-step. Therefore, the extended Lagrangian for charge equilibration is comprised of two auxiliary variables for the intermediate charges which serve as an initial guess for the real charges. We show that the iEL-SCF is able to reduce the number of SCF cycles by 50-80% of the original conjugate gradient self-consistent field solver as tested across diverse systems including water, ferric hydroxide, nitramine RDX, and hexanitrostilbene.

Development of an Advanced Force Field for Water Using Variational Energy Decomposition Analysis.

Given the piecewise approach to modeling intermolecular interactions for force fields, they can be difficult to parametrize since they are fit to data like total energies that only indirectly connect to their separable functional forms. Furthermore, by neglecting certain types of molecular interactions such as charge penetration and charge transfer, most classical force fields must rely on, but do not always demonstrate, how cancellation of errors occurs among the remaining molecular interactions accounted for such as exchange repulsion, electrostatics, and polarization. In this work we present the first generation of the (many-body) MB-UCB force field that explicitly accounts for the decomposed molecular interactions commensurate with a variational energy decomposition analysis, including charge transfer, with force field design choices that reduce the computational expense of the MB-UCB potential while remaining accurate. We optimize parameters using only a single water molecule and water cluster data up through pentamers, with no fitting to condensed phase data, and we demonstrate that high accuracy is maintained when the force field is subsequently validated against conformational energies of larger water cluster data sets, radial distribution functions of the liquid phase, and the temperature dependence of thermodynamic and transport water properties. We conclude that MB-UCB is comparable in performance to MB-Pol but is less expensive and more transferable by eliminating the need to represent short-ranged interactions through large parameter fits to high order polynomials.

Quantifying the two-state facilitated diffusion model of protein-DNA interactions.

The current report extends the facilitated diffusion model to account for conflict between the search and recognition binding modes adopted by DNA-binding proteins (DBPs) as they search DNA and subsequently recognize and bind to their specific binding site. The speed of the search dynamics is governed by the energetic ruggedness of the protein-DNA landscape, whereas the rate for the recognition process is mostly dictated by the free energy barrier for the transition between the DBP's search and recognition binding modes. We show that these two modes are negatively coupled, such that fast 1D sliding and rapid target site recognition probabilities are unlikely to coexist. Thus, a tradeoff occurs between optimizing the timescales for finding and binding the target site. We find that these two kinetic properties can be balanced to produce a fast timescale for the total target search and recognition process by optimizing frustration. Quantification of the facilitated diffusion model by including a frustration term enables it to explain several experimental observations concerning search and recognition speeds. The extended model captures experimental estimate of the energetic ruggedness of the protein-DNA landscape and predicts how various molecular properties of protein-DNA binding affect recognition kinetics. Particularly, point mutations may change the frustration and so affect protein association with DNA, thus providing a means to modulate protein-DNA affinity by manipulating the protein's association or dissociation reactions.

Smallest Archimedean Screw: Facet Dynamics and Friction in Multiwalled Nanotubes.

We identify a new material phenomenon, where minute mechanical manipulations induce pronounced global structural reconfigurations in faceted multiwalled nanotubes. This behavior has strong implications on the tribological properties of these systems and may be the key to understand the enhanced interwall friction recently measured for boron-nitride nanotubes with respect to their carbon counterparts. Notably, the fast rotation of helical facets in these systems upon coaxial sliding may serve as a nanoscale Archimedean screw for directional transport of physisorbed molecules.

Multiwalled nanotube faceting unravelled.

Nanotubes show great promise for miniaturizing advanced technologies. Their exceptional physical properties are intimately related to their morphological and crystal structure. Circumferential faceting of multiwalled nanotubes reinforces their mechanical strength and alters their tribological and electronic properties. Here, the nature of this important phenomenon is fully rationalized in terms of interlayer registry patterns. Regardless of the nanotube identity (that is, diameter, chirality, chemical composition), faceting requires the matching of the chiral angles of adjacent layers. Above a critical diameter that corresponds well with experimental results, achiral multiwalled nanotubes display evenly spaced extended axial facets whose number equals the interlayer difference in circumferential unit cells. Elongated helical facets, commonly observed in experiment, appear in nanotubes that exhibit small interlayer chiral angle mismatch. When the wall chiralities are uncorrelated, faceting is suppressed and outer layer corrugation, which is induced by the Moiré superlattice, is obtained in agreement with experiments. Finally, we offer an explanation for the higher incidence of faceting in multiwalled boron nitride nanotubes with respect to their carbon-based counterparts.

Coherent commensurate electronic states at the interface between misoriented graphene layers.

Graphene and layered materials in general exhibit rich physics and application potential owing to their exceptional electronic properties, which arise from the intricate π-orbital coupling and the symmetry breaking in twisted bilayer systems. Here, we report room-temperature experiments to study electrical transport across a bilayer graphene interface with a well-defined rotation angle between the layers that is controllable in situ. This twisted interface is artificially created in mesoscopic pillars made of highly oriented pyrolytic graphite by mechanical actuation. The overall measured angular dependence of the conductivity is consistent with a phonon-assisted transport mechanism that preserves the electron momentum of conduction electrons passing the interface. The most intriguing observations are sharp conductivity peaks at interlayer rotation angles of 21.8° and 38.2°. These angles correspond to a commensurate crystalline superstructure leading to a coherent two-dimensional (2D) electronic interface state. Such states, predicted by theory, form the basis for a new class of 2D weakly coupled bilayer systems with hitherto unexplored properties and applications.

Interlayer Potential for Graphene/h-BN Heterostructures.

We present a new force-field potential that describes the interlayer interactions in heterojunctions based on graphene and hexagonal boron nitride (h-BN). The potential consists of a long-range attractive term and a short-range anisotropic repulsive term. Its parameters are calibrated against reference binding and sliding energy profiles for a set of finite dimer systems and the periodic graphene/h-BN bilayer, obtained from density functional theory using a screened-exchange hybrid functional augmented by a many-body dispersion treatment of long-range correlation. Transferability of the parametrization is demonstrated by considering the binding energy of bulk graphene/h-BN alternating stacks. Benchmark calculations for the superlattice formed when relaxing the supported periodic heterogeneous bilayer provide good agreement with both experimental results and previous computational studies. For a free-standing bilayer we predict a highly corrugated relaxed structure. This, in turn, is expected to strongly alter the physical properties of the underlying monolayers. Our results demonstrate the potential of the developed force-field to model the structural, mechanical, tribological, and dynamic properties of layered heterostructures based on graphene and h-BN.

Inter-layer potential for hexagonal boron nitride.

A new interlayer force-field for layered hexagonal boron nitride (h-BN) based structures is presented. The force-field contains three terms representing the interlayer attraction due to dispersive interactions, repulsion due to anisotropic overlaps of electron clouds, and monopolar electrostatic interactions. With appropriate parameterization, the potential is able to simultaneously capture well the binding and lateral sliding energies of planar h-BN based dimer systems as well as the interlayer telescoping and rotation of double walled boron-nitride nanotubes of different crystallographic orientations. The new potential thus allows for the accurate and efficient modeling and simulation of large-scale h-BN based layered structures.

Robust Superlubricity in Graphene/h-BN Heterojunctions.

The sliding energy landscape of the heterogeneous graphene/h-BN interface is studied by means of the registry index. For a graphene flake sliding on top of h-BN, the anisotropy of the sliding energy corrugation with respect to the misfit angle between the two naturally mismatched lattices is found to reduce with the flake size. For sufficiently large flakes, the sliding energy corrugation is expected to be at least an order of magnitude lower than that obtained for matching lattices regardless of the relative interlayer orientation. Therefore, in contrast to the case of the homogeneous graphene interface where flake reorientations are known to eliminate superlubricty, here, a stable low-friction state is expected to occur. Our results mark heterogeneous layered interfaces as promising candidates for dry lubrication purposes.

Ultrahigh torsional stiffness and strength of boron nitride nanotubes.

We report the experimental and theoretical study of boron nitride nanotube (BNNT) torsional mechanics. We show that BNNTs exhibit a much stronger mechanical interlayer coupling than carbon nanotubes (CNTs). This feature makes BNNTs up to 1 order of magnitude stiffer and stronger than CNTs. We attribute this interlayer locking to the faceted nature of BNNTs, arising from the polarity of the B-N bond. This property makes BNNTs superior candidates to replace CNTs in nanoelectromechanical systems (NEMS), fibers, and nanocomposites.