Submicrosecond-Scale MD Simulations

• Submicrosecond-scale MD simulations are atomistic models that capture molecular events over hundreds of nanoseconds to microseconds, enabling observation of rare and long-timescale kinetics.
• Advances such as wafer-scale hardware, neural network potentials, and hybrid NNP/MM methods dramatically increase simulation speed and accuracy, overcoming classical integration bottlenecks.
• Enhanced sampling techniques and generative surrogate models further optimize computational cost and reliability, paving the way for routine submicrosecond-to-millisecond atomistic studies.

Submicrosecond-scale molecular dynamics (MD) simulations refer to the direct atomistic simulation of molecular systems over timescales reaching hundreds of nanoseconds to microseconds per trajectory. Achieving such durations is essential for capturing rare events, conformational changes, and long-timescale kinetics in biomolecular, chemical, and materials processes. However, realizing these timescales with sufficient accuracy and physical realism has presented a persistent computational challenge, addressed through advances in hardware, algorithmic sampling, machine learning potentials, and generative surrogates.

1. Historical Perspective and Computational Bottlenecks

The integration timestep in classical MD is restricted by the fastest molecular motions, typically bond vibrations, necessitating timesteps on the order of 1 femtosecond (fs). As a result, direct simulation of microsecond or longer trajectories can require over $10^{9}$ time steps, placing formidable demands on computational resources. For decades, improvements in hardware throughput yielded proportional gains in accessible simulation time, but this strong scaling has plateaued. Modern scale-out CPU and GPU architectures enhance spatial system size but not temporal coverage, due to communication overheads and the fixed sequential nature of explicit MD integration (Perez et al., 15 Nov 2024, Santos et al., 13 May 2024). Bespoke MD hardware solutions, such as Anton, allowed specialized acceleration but offered limited programmability and community impact.

2. Hardware Advances and Exascale Architectures

The recent introduction of wafer-scale platforms, exemplified by the Cerebras Wafer Scale Engine (WSE), has redefined strong scaling for atomistic MD. The WSE-2 integrates over 850,000 programmable cores in a monolithic silicon wafer, supporting up to 1.144 million MD steps per second for 200,000 atoms with Embedded Atom Method potentials; this performance is three orders of magnitude higher in simulated time per day than exascale GPUs (Perez et al., 15 Nov 2024, Santos et al., 13 May 2024). Each atom maps to an individual processor core, eliminating distributed-memory communication bottlenecks. Local memory and mesh-based message-passing enable atom-core and spatial locality preservation, while computation is load-balanced dynamically to reflect system evolution. These advances unlock direct simulation of millisecond processes and processes governed by slow kinetics (e.g., microstructure transformation in materials), previously unapproachable on general-purpose hardware.

<table> <tr><th>Platform</th><th>Atoms</th><th>Steps/sec</th><th>ms/day</th></tr> <tr><td>Cerebras WSE-2</td><td\>200,000</td><td\>1,144,000</td><td>~0.1</td></tr> <tr><td>Anton-3</td><td\>24,000</td><td\>980,000</td><td>~0.08</td></tr> <tr><td>Frontier (GPU)</td><td\>800,000</td><td\>1,530</td><td\>0.00009</td></tr> </table>

For neural network molecular dynamics (NNMD) with ab initio accuracy, innovations in intra-node MPI communication, threading, and kernel optimization (e.g., SVE-512 intrinsics, mixed precision GEMM) have pushed DeePMD-kit to 149 ns/day on 12,000 Fugaku nodes—over 30× faster than prior implementations, extending feasible simulation timescales into the millisecond regime (Li et al., 30 Oct 2024).

3. Hybrid, ML-Driven, and Data-Efficient Potentials

The use of machine-learned potentials—typically neural network potentials (NNPs)—has propelled accuracy toward chemical precision but at significant computational cost. Submicrosecond MD leveraging these models requires both algorithmic and implementation advances.

Hybrid Methods

The NNP/MM hybrid method partitions the simulation into an NNP region (e.g., ligand) and a molecular mechanics (MM) region (protein, solvent), inspired by QM/MM approaches but with greater efficiency. GPU-resident computational kernels and ensemble-parallelized NNP evaluations (e.g., TorchANI-2x with NNPOps) yield 5× speedups, allowing unbiased 1 μs trajectories in protein–ligand systems (Galvelis et al., 2022). Performance benchmarks confirm up to 74 ns/day for complex systems, with the bottleneck remaining in NNP inference. Limitations include element/charge coverage (e.g., only neutral ligands for ANI-2x) and the need for reduced time steps (2 fs) for stability.

Transfer Learning for High-Level PES

Transfer learning enables building CCSD(T)-quality potential energy surfaces (PES) for gas-phase MD by fine-tuning a low-level (e.g., HF or DFT) neural network PES with a small number (~100) of high-level ab initio points (Käser et al., 2023). This approach yields submicrosecond stable simulations for malonaldehyde with mean absolute errors of <0.06 kcal/mol (barriers) and <3 cm⁻¹ (frequencies). Simulation of ~1 μs (multiple 250 ps NVE trajectories) is feasible with <1% of the computational effort of direct ab initio MD.

Multiple Time-Step (MTS) and Distillation Strategies

For large systems and foundation NNPs, a dual-level MTS framework combines a short-range, fast distilled model (3.5 Å cutoff) for bonded interactions with a slow, full NNP evaluated every 3–6 fs in a RESPA-like integration loop. This reduces computational expense by 2.3–4× with negligible loss of accuracy for properties including free energies and protein structures (Cattin et al., 8 Oct 2025). Hydrogen mass repartitioning (HMR) supports even larger outer time steps with robust energy conservation.

4. Enhanced Sampling and Population-Based Algorithms

Enhanced sampling methods accelerate convergence without requiring direct propagation by small time steps throughout.

Population Annealing MD (PAMD)

PAMD combines parallel ensembles with resampling at stepwise cooling, borrowing from Monte Carlo population annealing. Stochastic thermostats decorrelate replicas after Boltzmann re-weighting, and physical observables are obtained as population averages. Unlike parallel tempering, PAMD scales linearly to thousands of replicas, utilizing all available hardware (Christiansen et al., 2018). For peptide folding, effective sampling equivalent to long continuous trajectories is achieved in significantly reduced wall-clock time.

Randomized Neighbor Schemes

The random-batch list (RBL) algorithm replaces full Verlet neighbor lists with two-level (core/shell) lists and stochastic batch estimation for shell interactions (Liang et al., 2021). Computational cost and storage fall by up to 10× without measurable loss in key observables, supporting longer simulations for dense or large systems by relieving the bottleneck of pair list construction.

5. Generative Surrogates and Latent Space Simulation

Several recent approaches move beyond explicit time integration, recasting the simulation as sequential or MCMC sampling in reduced or transformed spaces.

Markov State Models (MSMs) and synMD

Fine-grained MSMs, built with non-standard stratified clustering, enable efficient stochastic trajectory propagation on nanosecond–microsecond lag times, with backmapping to atomistic structures. This yields orders-of-magnitude acceleration over direct MD, with preservation of mean first passage times (MFPTs) and equilibrium populations up to the temporal resolution of the underlying MSM (Russo et al., 2022).

Latent Space Simulators (LSS) and Equation-Free Frameworks

Three-network architectures—(1) encoder for slow collective variables (SRV), (2) propagator (MDN), (3) all-atom decoder (cWGAN)—allow efficient simulation and reconstruction of long continuous trajectories (Sidky et al., 2020). Kinetics and thermodynamics match reference MD within statistical uncertainties, and sampling efficiency improves by ~10⁶×, as demonstrated for Trp-cage folding.

The Learning Effective Dynamics (LED) framework combines MDN-autoencoding with LSTM-MDN non-Markovian time propagation, allowing sampling of long timescale phenomena at three orders of magnitude lower cost per simulated time than MD, without loss of accuracy in thermodynamic distributions or transition rates (Vlachas et al., 2021).

6. Probabilistic and Diffusion-Based Generative Models

Normalizing flow models, as in Timewarp, are trained to generate large time-step proposals (e.g., 100 ps, $10^{5-6}$ fs) in MCMC schemes targeting the Boltzmann distribution. These models generalize to new small peptides and achieve up to 600× acceleration in mapping metastable states compared to MD (Klein et al., 2023). Score Dynamics (SD), based on graph neural networks and diffusion models, samples updates over 10 ps intervals—four orders of magnitude larger than explicit integrators—for small solutions, preserving equilibrium and kinetic observables with 80–180× wall-clock speedup (Hsu et al., 2023). Both frameworks are currently constrained by training data requirements and performance in memory-dominated or topologically complex systems.

7. Future Directions and Ongoing Limitations

Current submicrosecond-capable methods each face domain-specific trade-offs:

• Hybrid NNP/MM approaches are limited by physics coverage of NNPs (element, charge) and ensemble evaluation cost, but software and model improvements are ongoing (Galvelis et al., 2022).
• Generative and latent variable models are bounded by the diversity and quality of training data and often cannot extrapolate to unseen thermodynamic states or topologies without retraining (Sidky et al., 2020, Hsu et al., 2023).
• Score-based diffusion models require substantial training data; stability and conservation over very long timescales present open challenges (Hsu et al., 2023).
• Wafer-scale and node-based architectures for direct MD are general but currently optimized primarily for structured potentials (e.g., EAM or DeePMD) and not yet coupled to the full range of ML-driven ab initio models (Perez et al., 15 Nov 2024, Li et al., 30 Oct 2024, Santos et al., 13 May 2024).

A plausible implication is that as these platforms and hybrid algorithms mature—combining wafer-scale hardware, efficient ML potentials, and advanced sampling—they may close the gap between simulation throughput, accuracy, and accessible process timescales, making direct submicrosecond–millisecond atomistic MD routine across disciplines.