Small ionic radii limit time step in Martini 3 molecular dynamics simulations

TL;DR

This paper investigates how small ionic bead sizes in Martini 3 force field limit simulation time steps, and proposes methods to increase time steps without compromising accuracy, enhancing simulation efficiency.

Contribution

It introduces a kinetic model of MD instabilities related to ion representation and demonstrates how increasing ionic mass or size allows longer time steps in Martini 3 simulations.

Findings

Time step limited to 25 fs with tiny ion beads in Martini 3.

Increasing ionic mass or bead size enables up to 40 fs time steps.

Ion-water interactions are primary source of simulation instability.

Abstract

Among other improvements, the Martini 3 coarse-grained force field provides a more accurate description of the solvation of protein pockets and channels through the consistent use of various bead types and sizes. Here, we show that the representation of Na$^+$ and Cl$^-$ ions as "tiny" (TQ5) beads limits the accessible time step to 25 fs. By contrast, with Martini 2, time steps of 30-40 fs were possible for lipid bilayer systems without proteins. This limitation is relevant for, e.g., phase separating lipid mixtures that require long equilibration times. We derive a quantitative kinetic model of time-integration instabilities in molecular dynamics (MD) as a function of time step, ion concentration and mass, system size, and simulation time. With this model, we demonstrate that ion-water interactions are the main source of instability at physiological conditions, followed closely by ion-ion interactions. We show that increasing the ionic masses makes it possible to use time steps up to 40 fs with minimal impact on static equilibrium properties and on dynamical quantities such as lipid and ion diffusion coefficients. Increasing the size of the bead representing the ions (and thus changing their hydration) also permits longer time steps. The use of larger time steps in Martini 3 simulations results in a more efficient exploration of configuration space. The kinetic model of MD simulation crashes can be used to determine the maximum allowed time step whenever sampling efficiency is critical.