Calculation of the electrostatic potential of lipid bilayers from molecular dynamics simulations: methodological...

The electrostatic properties of lipid membranes are of profound importance as they are directly associated with membrane potential and, consequently, with numerous membrane-mediated biological phenomena. Here we address a number of methodological issues related to the computation of the electrostatic potential from atomic-scale molecular dynamics simulations of lipid bilayers. We discuss two slightly different forms of Poisson equation that are normally used to calculate the membrane potential: (i) a classical form when the potential and the electric field are chosen to be zero on one of the sides of a simulation box and (ii) an alternative form, when the potential is set to be the same on the opposite sides of a simulation box. Both forms differ by a position-dependent correction term, which has been shown to be proportional to the overall dipole moment of a bilayer system (for neutral systems). For symmetric bilayers we demonstrate that both approaches give essentially the same potential profiles, provided that simulations are long enough (a production run of at least 100 ns is required) and that fluctuations of the center of mass of a bilayer are properly accounted for. In contrast, for asymmetric lipid bilayers, the second approach is no longer appropriate due to a nonzero net dipole moment across a simulation box with a single asymmetric bilayer. We demonstrate that in this case the electrostatic potential can adequately be described by the classical form of Poisson equation, provided that it is employed in conjunction with tin-foil boundary conditions, which exactly balance a nonzero surface charge of a periodically replicated multibilayer system. Furthermore, we show that vacuum boundary conditions give qualitatively similar potential profiles for asymmetric lipid bilayers as compared to the conventional periodic boundaries, but accurate determination of the transmembrane potential difference is then hindered due to detachment of some water dipoles from bulk aqueous solution to vacuum.