fix pimd/langevin: An Efficient Implementation of Path Integral Molecular Dynamics in LAMMPS

TL;DR

This work introduces fix pimd/langevin, an efficient, MPI-parallel PIMD implementation embedded in LAMMPS, enabling large-scale simulations of nuclear quantum effects with modern ML interatomic potentials. Grounded in the normal-mode PIMD framework, it supports NVE/NVT/NPT ensembles via PILE_L and BZP, two-level MPI parallelization, and optional bosonic PIMD, and it seamlessly handles periodic boundaries. The authors validate correctness against i-PI using liquid water and demonstrate several-fold speedups, coupled with strong and weak scaling that reveal robust performance on massively parallel hardware, including DP-driven simulations with long-range electrostatics (DPLR). They also provide practical usage guidance, including long-range electrostatics handling, example workflows, and a comprehensive feature set, underscoring the practical impact for high-accuracy, large-system quantum simulations in materials and chemistry. Overall, fix pimd/langevin represents a significant step toward efficient, scalable PIMD integrated directly into a widely-used MD engine, enabling ab initio-like quantum effects in large MLIP-enabled systems.

Abstract

Path integral molecular dynamics (PIMD), which maps a quantum particle onto a fictitious classical system of ring polymers and propagates the "beads" of this extended classical system using molecular dynamics, is widely used to capture nuclear quantum effects (NQEs) in molecular simulations. Accurate PIMD calculations typically require a large number of beads and are therefore computationally demanding. While software packages such as i-PI offer comprehensive PIMD functionality, the high efficiency of simulations driven by machine learning interatomic potentials, such as Deep Potential (DP), calls for more efficient PIMD implementations that fully exploit modern massively parallel supercomputers. Here we present fix pimd/langevin, an efficient PIMD implementation in LAMMPS that supports commonly used features and leverages the Message Passing Interface architecture of LAMMPS to achieve high computational efficiency. We demonstrate the usage and validate the correctness of our code using liquid water as a representative example, and provide a comprehensive overview of the supported features. Then we discuss several important technical aspects of the implementation. Using DP simulations of water as a benchmark, we show that our implementation achieves several-fold acceleration compared to i-PI. Finally, we report strong and weak scaling results that demonstrate the favorable parallel performance of our code.

Figures (11)

• Figure 1: Centroid-virial quantum kinetic energy estimator $K_{\mathrm{CV}}$. Left: trajectory of the quantity (solid line) and its cumulative average (dashed line). Right: corresponding distribution. Blue lines correspond to LAMMPS results and red lines to i-PI results.
• Figure 2: Temperature $T$ of the classical ring polymer system. The temperature printed by LAMMPS's temp should be divided by the number of beads for a proper comparison. Left: trajectory of the quantity (solid line) and its cumulative average (dashed line). Right: corresponding distribution. Blue lines correspond to LAMMPS results and red lines to i-PI results.
• Figure 3: Centroid-virial pressure estimator $P_{\mathrm{CV}}$. Left: trajectory of the quantity (solid line) and its cumulative average (dashed line). Right: corresponding distribution. Blue lines correspond to LAMMPS results and red lines to i-PI results.
• Figure 4: The density $\rho$ of the system. Left: trajectory of the quantity (solid line) and its cumulative average (dashed line). Right: corresponding distribution. Blue lines correspond to LAMMPS results and red lines to i-PI results.
• Figure 5: The radial distribution functions (RDFs) of liquid water from PIMD simulations for 128 H$_2$O molecules at 300 K and 1 bar. The blue line corresponds to LAMMPS results and the red line to i-PI results.