First principles-based multiscale atomistic methods for input into first principles nonequilibrium transport across interfaces

Abstract

This issue of PNAS features "nonequilibrium transport and mixing across interfaces," with several papers describing the nonequilibrium coupling of transport at interfaces, including mesoscopic and macroscopic dynamics in fluids, plasma, and other materials over scales from microscale to celestial. Most such descriptions describe the materials in terms of the density and equations of state rather than specific atomic structures and chemical processes. It is at interfacial boundaries where such atomistic information is most relevant. However, there is not yet a practical way to couple these phenomena with the atomistic description of chemistry. The starting point for including such information is the quantum mechanics (QM). However, practical QM calculations are limited to a hundred atoms for dozens of picoseconds, far from the scales required to inform the continuum level with the proper atomistic description. To bridge this enormous gap, we need to develop practical methods to extend the scale of the atomistic simulation by several orders of magnitude while retaining the level of QM accuracy in describing the chemical process. These developments would enable continuum modeling of turbulent transport at interfaces to incorporate the relevant chemistry. In this perspective, we will focus on recent progress in accomplishing these extensions in first principles-based atomistic simulations and the strategies being pursued to increase the accuracy of very large scales while dramatically decreasing the computational effort.

Keywords: electron force field; molecular dynamics; multiscale simulation; quantum mechanics; reactive force fields.

Fig. 1.

Comparison of ReaxFF with QM (on the level of B3LYP/6-311G*) for decomposition of RDX molecule via all four reaction pathways (3).

Fig. 2.

QM polarization as a dipole is brought up to cyclohexane compared with PQEq, QEq, and common FFs, such as optimized potential for liquid simulations (OPLS). We see that PQEq and QEq agree well with QM, while static charge models do not.

Fig. 3.

(A) Comparison of equation of state for H₂ crystal from QM and experiment. (B) Inset shows QM-derived nonbond interactions for He, Ne, Ar, Kr, Xe, and Rn crystals. This is a demonstration that they scale to a universal NB curve (black solid curve).

Fig. 4.

Illustration of the O-H bond in water, which shows how to define the RexPoN BE that, added to the QM-based NB, exactly fits highest-level QM. NBrep, the repulsion part of nonbonded interaction; NBatt, the attractive part of nonbonded interaction.

Fig. 5.

Application of eFF to Auger decay-induced etching of the C₁₉₇H₁₁₂ diamondoid. (A) Distance of valence electrons from the core hole. The figure shows the green electron filling the core hole, the red electron being ejected (and trapped after 20 femtoseconds; not shown), and the blue and purple electrons being excited. Inset depicts color-coded trajectories of valence electrons during a core-hole filling event. (B) Charges-associated atomic composition of fragments separated from the nanoparticle (NP) after 50 femtoseconds.

Fig. 6.

The half-life time for H₂O product formation predicted from 798 K to 2,698 K using aARRDyn. The time for forming one-half the H₂O products at 798 K is 538 seconds (just 1.3 million time steps).

Fig. 7.

Simulation of brittle fracture in B4C using the ReaxFF trained with QM. (A) The periodic cell has 200,000 atoms. (B) Stress–strain relationship as the system is sheared to brittle failure. (C) The successive processes of twinning, amorphous band formation, cavitation, and crack formation: (B) twinning initiation, (C) twins grow, (D) amorphous band initiates from twins, (E) cavitation within amorphous band, and (F) cavity grows into a crack and then massive failure.

Fig. 8.

Shock wave passing through a nonplanar sawtooth interface of PBX. (A) PBX model. (B) One-dimensional temperature profile at 6.5 picoseconds showing that the hot spot is suppressed by decreasing the density of the polymer. (C) Two-dimensional temperature profile showing hot spot formation at the asperity, which remains in place as the shock proceeds.

Fig. 9.

(A) ReaxFF-simulated TiO₂ crater radius versus ice cluster impactor radius is consistent with experimental data from Kearsley et al. (56) associated with NASA’s Stardust mission. Atomistic crater profiles are shown as Insets. The ratio of crater radius (r_c) to cluster radius (r_p) is linear over the range of impact velocities and its slope (s) provides an indication of increasing or decreasing damage [i.e., at encounter (E)5, s 1⁄4 3:6; therefore, larger-diameter damage is expected than at E3 with s 1⁄4 3:7]. Kearsley et al. (56) report s 1⁄4 2:2–4:63. (B) Sublimated titanium atoms (in orange) from INMS walls would potentially affect its accuracy in reading oxygen-containing species.

Fig. 10.

Observed (solid) and analytically reproduced (dashed) number of species as a function of time for microcanonical ensemble (NVE) RMD simulation equilibrated at (A) 1,000 K and (B) 2,000 K.