Lipid bilayers driven to a wrong lane in molecular dynamics simulations by subtle changes in long-range electrostatic interactions

We provide compelling evidence that different treatments of electrostatic interactions in molecular dynamics simulations may dramatically affect dynamic properties of lipid bilayers. To this end, we consider a fully hydrated pure dipalmitoylphosphatidylcholine bilayer through 50-ns molecular dynamics simulations and study various dynamic properties of individual lipids in a membrane, including the velocity autocorrelation function, the lateral and rotational diffusion coefficients, and the autocorrelation function for the area per molecule. We compare the results based on the Particle-Mesh Ewald (PME) and reaction field (RF) techniques with those obtained by an approach where the electrostatic interactions are truncated at r(cut) = 1.8, 2.0, and 2.5 nm. We find that the overall performance of PME is very good; its results are consistent with the expected behavior. The RF method performs rather well, too, despite certain inherent problems and the fact that its results differ from those obtained by PME. Nevertheless, the largest differences are found for the truncation methods, for which all examined truncation methods lead to results distinctly different from those obtained by PME. The lateral diffusion coefficients obtained by PME and truncation at 1.8 nm, for example, differ by a factor of 10, while the PME results are consistent with experimental values. The observed deviations can be interpreted in terms of artificial ordering due to truncation and highlight the important role of electrostatic interactions in the dynamics of systems composed of lipids and other biologically relevant molecules such as proteins and DNA.