A generalized reaction field method for molecular dynamics simulations

TLDR

Molecular dynamics simulations of ionic systems require inclusion of long‑range electrostatic forces. The authors propose an analytically derived expression for long‑range electrostatic forces based on the Poisson–Boltzmann equation outside a spherical cutoff, readily implementable in simulation codes. They solve the linearized Poisson–Boltzmann equation in a spherical region, derive potential and force expressions at the sphere center, implement these as a generalized reaction field in MD, and test the method on a periodic 2127‑molecule SPC water box with 40 NaCl ions, comparing results to simple cutoff, cutoff plus PB reaction field, and Ewald summation. Using the Poisson–Boltzmann generalized reaction field markedly improves solution structure over simple cutoff, without additional computational cost.

Abstract

Molecular dynamics simulations of ionic systems require the inclusion of long-range electrostatic forces. We propose an expression for the long-range electrostatic forces based on an analytical solution of the Poisson–Boltzmann equation outside a spherical cutoff, which can easily be implemented in molecular simulation programs. An analytical solution of the linearized Poisson–Boltzmann (PB) equation valid in a spherical region is obtained. From this general solution special expressions are derived for evaluating the electrostatic potential and its derivative at the origin of the sphere. These expressions have been implemented for molecular dynamics (MD) simulations, such that the surface of the cutoff sphere around a charged particle is identified with the spherical boundary of the Poisson–Boltzmann problem. The analytical solution of the Poisson–Boltzmann equation is valid for the cutoff sphere and can be used for calculating the reaction field forces on the central charge, assuming a uniform continuum of given ionic strength beyond the cutoff. MD simulations are performed for a periodic system consisting of 2127 SPC water molecules with 40 NaCl ions (1 molar). We compare the structural and dynamical results obtained from MD simulations in which the long range electrostatic interactions are treated differently; using a cutoff radius, using a cutoff radius and a Poisson–Boltzmann generalized reaction field force, and using the Ewald summation. Application of the Poisson–Boltzmann generalized reaction field gives a dramatic improvement of the structure of the solution compared to a simple cutoff treatment, at no extra computational cost.