next up previous contents
Next: Analysis of Simulation Data Up: Practicalities Previous: Periodic Boundary Conditions   Contents


Long Range Interactions

Non-bonded interactions can be divided into two classes; short and long range interactions. Formally a force is defined to be short ranged if it decreases with distance quicker than $r^{-d}$ where $d$ is the dimensionality of the system (usually 3). Short ranged interactions are commonly dealt with by imposing a cut-off to the potential $V(r)$, $r_c$, beyond which $V(r)$ is set to 0. This can be justified by the speed at which a short range force decays with distance. As long as $r_c$ is chosen to be sufficiently large, the cut-off will only imposes a slight perturbation on the system (although $r_c$ has to be less than half the box length to be consistent with the minimum image convention). Correction formulae can be applied to quantities such as the pressure or energy to correct for this. Long range forces present more of a problem. As they are infinite ranged, the simulation cell and all its periodic images must be considered. In principle for a large enough system screening by neighbouring molecules would diminish the effect of the potential. However, this would occur over a range of several tens or hundreds of nanometres while current simulations typically have box lengths of the order of nanometres. Thus to account for long range interactions in systems that can be reasonably simulated we must consider the effect of the periodic images. In this case the Coulombic interaction for a set of point charges is
\begin{displaymath}
U_{coul}=\frac{1}{4\pi \epsilon_0}\sum_{\mathbf{n}}
\l...
... q_i q_j \vert \mathbf{r}_{ij}+\mathbf{n}\vert^{-1}
\right)
\end{displaymath} (3.46)

where $\mathbf{n}$ are the lattice vectors $\mathbf{n}=(n_x L_x, n_y L_y,
n_zL_z)$. This sum is conditionally convergent; the result depends on the order the terms are summed in. The Ewald sum [127,115] is commonly used to model the effect of long range forces. It was originally developed for the study of ionic crystals. The Ewald sum decomposes the sum in (3.46) into two rapidly convergent sums
\begin{displaymath}
\sum\frac{1}{r}=\sum_m \frac{F(m)}{r}+\sum_n \frac{1-F(n)}{r}.
\end{displaymath} (3.47)

Physically the Ewald sum works by surrounding each point charge in the system by a charge distribution of equal magnitude and opposite sign. This is commonly taken to be a Gaussian distribution, although this choice is arbitrary. The counter-charge screens the original potential making it short ranged. This is then summed in real space. Then a second imaginary charge distribution of opposite sign to the first (and of the same sign as the point charges) is added to cancel out the screening charge. As this screening distribution is a smooth function, its Fourier transform is rapidly convergent. Thus this second is summed in reciprocal space. Two further terms are present in the Ewald sum. The first is the self term. This cancels the interaction between a point charge and its own screening distribution and is a constant. The second is the surface term.This accounts for the dipolar layer that appears at the surface of a sphere in a vacuum. The final potential energy is then
$\displaystyle U_{coul}$ $\textstyle =$ $\displaystyle \frac{1}{2}\sum_{i=1}^{N}\sum_{j=1}^{N}q_iq_j\frac{\mbox{erfc}
(\alpha\vert
\mathbf{r}_{ij}+\mathbf{n}\vert)}{\vert\mathbf{r}_{ij}+\mathbf{n}\vert}$  
  $\textstyle +$ $\displaystyle \frac{1}{\pi L^3}\sum_{\mathbf{k}\not=\mathbf{0}} q_i q_j
\frac{4...
...2}{k^2}\exp\left(-\frac{k^2}{4\alpha^2}\right)
\cos(\mathbf{k}.\mathbf{r}_{ij})$  
  $\textstyle -$ $\displaystyle \frac{\alpha}{\pi^{\frac{1}{2}}}\sum_{i=1}^{N}q_i^2 +\frac{2\pi}
{3L^3}\left\vert \sum_{i=1}^{N}q_i \mathbf{r}_i\right\vert^2.$ (3.48)

The Ewald sum is possibly the most common method for evaluating long-range interactions in simulations. However, it can be expensive for large systems due to its $\mathcal{O}(N^{\frac{3}{2}})$ scaling. Thus several, cheaper, alternative methods have been proposed. One such technique is the Particle-Particle Particle-Mesh method (PPPM) [128] which scales as $\mathcal{O}(N\ln N)$. In this interactions at short range are directly calculated. Interactions at long range are evaluated by interpolating the charges onto a mesh. The Poisson equation can then be solved on this mesh to calculate the potential and forces. There are several variants of this technique that use different methods for discretizing the charges onto a mesh [129]. Another set of alternatives to the Ewald sum are the Fast Multipole Methods (FMM). These work by breaking the system into cells. The potential and force on each molecule are evaluated by considering the interaction between the molecule and each of these cells, where the distribution of charges in each of these cells are approximated by a multipole expansion. For a large enough system size this has been shown to scale as $\mathcal{O}(N)$, although due to a large prefactor, these methods only become favourable for systems of about 100000 particles [130]. Another alternative to the Ewald sum is the reaction field method [131,115]. Here contributions from molecules within a cavity of radius $r$ are calculated explicitly. Molecules outside this are considered to form a dielectric continuum producing a reaction field within the cavity.
next up previous contents
Next: Analysis of Simulation Data Up: Practicalities Previous: Periodic Boundary Conditions   Contents
Dr S J Clark
2003-01-30