Huai Sun

ReaxFF-lg: Correction of the ReaxFF Reactive Force Field for London Dispersion, with Applications to the Equations of State for Energetic Materials

The practical levels of density functional theory (DFT) for solids (LDA, PBE, PW91, B3LYP) are well-known not to account adequately for the London dispersion (van der Waals attraction) so important in molecular solids, leading to equilibrium volumes for molecular crystals ~10-15% too high. The ReaxFF reactive force field is based on fitting such DFT calculations and suffers from the same problem. In the paper we extend ReaxFF by adding a London dispersion term with a form such that it has low gradients (lg) at valence distances leaving the already optimized valence interactions intact but behaves as 1/R(6) for large distances. We derive here these lg corrections to ReaxFF based on the experimental crystal structure data for graphite, polyethylene (PE), carbon dioxide, and nitrogen and for energetic materials: hexahydro-1,3,5-trinitro-1,3,5-s-triazine (RDX), pentaerythritol tetranitrate (PETN), 1,3,5-triamino-2,4,6-trinitrobenzene (TATB), and nitromethane (NM). After this dispersion correction the average error of predicted equilibrium volumes decreases from 18.5 to 4.2% for the above systems. We find that the calculated crystal structures and equation of state with ReaxFF-lg are in good agreement with experimental results. In particular, we examined the phase transition between α-RDX and γ-RDX, finding that ReaxFF-lg leads to excellent agreement for both the pressure and volume of this transition occurring at ~4.8 GPa and ~2.18 g/cm(3) density from ReaxFF-lg vs 3.9 GPa and ~2.21 g/cm(3) from experiment. We expect ReaxFF-lg to improve the descriptions of the phase diagrams for other energetic materials.

Mechanism and Kinetics for the Initial Steps of Pyrolysis and Combustion of 1,6-Dicyclopropane-2,4-hexyne from ReaxFF Reactive Dynamics

We report the kinetic analysis and mechanism for the initial steps of pyrolysis and combustion of a new fuel material, 1,6-dicyclopropane-2,4-hexyne, that has enormous heats of pyrolysis and combustion, making it a potential high-energy fuel or fuel additive. These studies employ the ReaxFF force field for reactive dynamics (RD) simulations of both pyrolysis and combustion processes for both unimolecular and multimolecular systems. We find that both pyrolysis and combustion initiate from unimolecular reactions, with entropy-driven reactions being most important in both processes. Pyrolysis initiates with extrusion of an ethylene molecule from the fuel molecule and is followed quickly by isomerization of the fuel molecule, which induces additional radicals that accelerate the pyrolysis process. In the combustion process, we find three distinct mechanisms for the O(2) attack on the fuel molecule: (1) attack on the cyclopropane, ring expanding to form the cyclic peroxide which then decomposes; (2) attack onto the central single bond of the diyne which then fissions to form two C(5)H(5)O radicals; (3) attack on the alkyne-cyclopropane moiety to form a seven-membered ring peroxide which then decomposes. Each of these unimolecular combustion processes releases energy that induces additional radicals to accelerate the combustion process. Here oxygen has major effects both as the radical acceptor and as the radical producer. We extract both the effective activation energy and the effective pre-exponential factor by kinetic analysis of pyrolysis and combustion from these ReaxFF simulations. The low value of the derived effective activation energy (26.18 kcal/mol for pyrolysis and 16.40 kcal/mol for combustion) reveals the high activity of this fuel molecule.

COMPASS II: extended coverage for polymer and drug-like molecule databases

The COMPASS II force field has been developed by extending the coverage of the COMPASS force field (J Phys Chem B 102(38):7338–7364, 1998) to polymer and drug-like molecules found in popular databases. Using a fragmentation method to systematically construct small molecules that exhibit key functional groups found in these databases, parameters applicable to database compounds were efficiently obtained. Based on the same parameterization paradigm as used in the development of the COMPASS force field, new parameters were derived by a combination of fits to quantum mechanical data for valence parameters and experimental liquid and crystal data for nonbond parameters. To preserve the quality of the original COMPASS parameters, a quality assurance suite was used and updated to ensure that additional atom-types and parameters do not interfere with the existing ones. Validation against molecular properties, liquid and crystal densities, and enthalpies, demonstrates that the quality of COMPASS is preserved and the same quality of prediction is achieved for the additional coverage.

Development of a ReaxFF reactive force field for ettringite and study of its mechanical failure modes from reactive dynamics simulations.

Ettringite is a hexacalcium aluminate trisulfate hydrate mineral that forms during Portland cement hydration. Its presence plays an important role in controlling the setting rate of the highly reactive aluminate phases in cement paste and has also been associated with severe cracking in cured hardened cement. To understand how it forms and how its properties influence those of hardened cement and concrete, we have developed a first-principles-based ReaxFF reactive force field for Ca/Al/H/O/S. Here, we report on the development of this ReaxFF force field and on its validation and application using reactive molecular dynamics (RMD) simulations to characterize and understand the elastic, plastic, and failure response of ettringite at the atomic scale. The ReaxFF force field was validated by comparing the lattice parameters, pairwise distribution functions, and elastic constants of an ettringite crystal model obtained from RMD simulations with those from experiments. The predicted results are in close agreement with published experimental data. To characterize the atomistic failure modes of ettringite, we performed stress-strain simulations to find that Ca-O bonds are responsible for failure of the calcium sulfate and tricalcium aluminate (C3A) column in ettringite during uniaxial compression and tension and that hydrogen bond re-formation during compression induces an increase in plastic strain beyond the materials stress-strain proportionality limit. These results provide essential insight into understanding the mechanistic role of this mineral in cement and concrete degradation, and the ReaxFF potential developed in this work serves as a fundamental tool to further study the kinetics of hydration in cement and concrete.

Adaptive Accelerated ReaxFF Reactive Dynamics with Validation from Simulating Hydrogen Combustion

We develop here the methodology for dramatically accelerating the ReaxFF reactive force field based reactive molecular dynamics (RMD) simulations through use of the bond boost concept (BB), which we validate here for describing hydrogen combustion. The bond order, undercoordination, and overcoordination concepts of ReaxFF ensure that the BB correctly adapts to the instantaneous configurations in the reactive system to automatically identify the reactions appropriate to receive the bond boost. We refer to this as adaptive Accelerated ReaxFF Reactive Dynamics or aARRDyn. To validate the aARRDyn methodology, we determined the detailed sequence of reactions for hydrogen combustion with and without the BB. We validate that the kinetics and reaction mechanisms (that is the detailed sequences of reactive intermediates and their subsequent transformation to others) for H2 oxidation obtained from aARRDyn agrees well with the brute force reactive molecular dynamics (BF-RMD) at 2498 K. Using aARRDyn, we then extend our simulations to the whole range of combustion temperatures from ignition (798 K) to flame temperature (2998K), and demonstrate that, over this full temperature range, the reaction rates predicted by aARRDyn agree well with the BF-RMD values, extrapolated to lower temperatures. For the aARRDyn simulation at 798 K we find that the time period for half the H2 to form H2O product is ∼538 s, whereas the computational cost was just 1289 ps, a speed increase of ∼0.42 trillion (10(12)) over BF-RMD. In carrying out these RMD simulations we found that the ReaxFF-COH2008 version of the ReaxFF force field was not accurate for such intermediates as H3O. Consequently we reoptimized the fit to a quantum mechanics (QM) level, leading to the ReaxFF-OH2014 force field that was used in the simulations.

Urea: an ab initio and force field study of the gas and solid phases.

We have studied the gaseous and solid phases of urea using both quantum mechanics calculation and force field simulation methods. Our ab initio calculations confirmed experimental observations that urea structure is planar in the crystal, but nonplanar in the gas phase. Based on electron structure analysis, we suggest that the significant difference between these two structures in different environments can be qualitatively explained by two resonance structures. The planar structure is more polarized than the nonplanar one, and the former is stabilized in the solid phases due to strong electrostatic interactions. We found classical force field method is incapable to represent such strong polarization effect. Using molecular dynamics simulations with a force field optimized for condensed phases, we calculated the crystalline structures of urea in the temperature range of 12 to 293 K. The densities as well as cell parameters are within 2% deviation from the experimental data in the temperature range.

Classic Force Field for Predicting Surface Tension and Interfacial Properties of Sodium Dodecyl Sulfate

A classical force field capable for accurately predicting surface tensions, surface concentration, and other interfacial properties is reported for sodium dodecyl sulfate (SDS). This force field is proposed by combining parameters from well established force fields for components of the air/SDS/water interface and optimized by adjusting the van der Waals diameter of sulfate oxygen to fit experimental data of surface tension. The force field parameters are transferable; as good agreement with experimental data is achieved from independent calculations on activity derivatives of aqueous solution of sodium methyl sulfate using the Kirkwood-Buff theory. The adjusted parameter effectively modulates the electrostatic interactions of solvated ions in solutions. This modification has a strong impact to surface tension and location and mobility of sodium cations but minimal impact to properties such as density profiles of bulk phase and conformations and orientations of surfactant chain, for which consistent results compared with previous studies are obtained.

Stability of Newton Black Films Under Mechanical Stretch – A Molecular Dynamics Study

The stability of Newton black films (NBFs) under lateral mechanical stretch was investigated by nonequilibrium molecular dynamics (NEMD) simulations using force field parameters validated by accurate prediction of surface tensions. The applied strains accelerated film ruptures, enabling efficient measurements of the critical thicknesses of the films. Two representative surfactants, sodium dodecyl sulfate (SDS) for ionic surfactant and pentaethylene glycol monododecyl ether (C12EO5H) for nonionic surfactant, were investigated and compared. The predicted critical thickness of C12EO5H-coated film is smaller than that of the SDS-coated film, which is consistent with previously reported experimental observations. Our simulation results show that while the two surfactant-coated films exhibit similar dynamic properties attributed to the Marangoni-Gibbs effect, their surface structural characteristics are quite different. Consequently the two films demonstrate distinct rupture mechanisms in which rupture starts at uncovered water domains in the SDS-coated film, but at lateral surfactant/water interfaces in the C12EO5H-coated film. Our findings provide new insights into the stabilization mechanisms of NBFs and will facilitate the design and development of new films with improved properties.

Incorporating magnesium and calcium cations in porous organic frameworks for high-capacity hydrogen storage.

We propose incorporating a bi-functional group consisting of magnesium or calcium cations and a 1,2,4,5-benzenetetroxide anion (C6H2O4(4-)) in porous materials to enhance the hydrogen storage capacity. The C6H2O4M2 bifunctional group is highly stable and polarized, and each group provides 18 (M = Mg) or 22 (M = Ca) binding sites for hydrogen molecules with an average binding energy of ca. 10 kJ mol(-1) per hydrogen molecule based on RIMP2/ TZVPP calculations. Two porous materials (PAF-Mg or PAF-Ca) constructed with the bi-functional groups show remarkable improvement in hydrogen uptakes at normal ambient conditions. At 233 K and 10 MPa, the predicted gravimetric uptakes are 6.8 and 6.4 wt% for PAF-Mg and PAF-Ca respectively. This work reveals that fabricating materials with large numbers of binding sites and relatively low binding energies is a promising approach to achieve high capacity for on-board storage of hydrogen.

Hierarchical atom type definitions and extensible all-atom force fields

The extensibility of force field is a key to solve the missing parameter problem commonly found in force field applications. The extensibility of conventional force fields is traditionally managed in the parameterization procedure, which becomes impractical as the coverage of the force field increases above a threshold. A hierarchical atom‐type definition (HAD) scheme is proposed to make extensible atom type definitions, which ensures that the force field developed based on the definitions are extensible. To demonstrate how HAD works and to prepare a foundation for future developments, two general force fields based on AMBER and DFF functional forms are parameterized for common organic molecules. The force field parameters are derived from the same set of quantum mechanical data and experimental liquid data using an automated parameterization tool, and validated by calculating molecular and liquid properties. The hydration free energies are calculated successfully by introducing a polarization scaling factor to the dispersion term between the solvent and solute molecules.