Numerical Error Propagation in Contemporary Molecular Dynamics Simulations of Lithium Battery Components

Molecular dynamics (MD) simulations rely on accurate numerical integration of the equations of motion, where the choice of the timestep (Δt) critically affects stability and precision. Rather than aiming to recover an exact atomic trajectory, second-order integrators, such as those implemented in LAMMPS with the Nosé–Hoover thermostat, yield a global error scaling as O(Δt2). However, practical results deviate from this ideal behavior. In this work, we systematically analyze the effect of Δt on the accuracy of MD simulations using a polarizable force field applied to battery-relevant systems. Different errors were quantified for a range of timesteps. We evaluate whether timestep-induced numerical deviations affect physically meaningful observables, including diffusion coefficients and radial distribution functions, and examine the role of different integration schemes in solid and liquid components. The results show that decreasing Δt beyond the stability threshold leads to only marginal improvements in accuracy while significantly increasing computational cost at medium and long timescales. Conversely, excessively large Δt values produce numerical instability and unphysical behavior. These findings indicate that the expected ideal error scaling does not directly translate into practical accuracy gains and highlight the need for balanced timestep selection based on physical robustness rather than trajectory convergence alone.