Multilevel Summation Method for Calculating Electrostatics
Overview
Figure 0  Illustration of a semiperiodic membrane simulation using MSM. 
Significant longrange electrostatic interactions arise in many biomolecular systems, such as negatively charged DNA and RNA, polar or charged membranes, ion channels, and electrostatic steering of proteinprotein and enzymesubstrate association. Accordingly, electrostatic interactions need to be accurately represented in molecular modeling calculations.
The evaluation of electrostatic interactions can, through the use of the fast Fourier transform (FFT), be approximated with controlled accuracy by two finite range calculations, one in real space and one in reciprocal space. The particlemesh Ewald (PME) method implemented in NAMD is based on the FFT. However, the use of an FFTbased method, such as PME, has two drawbacks: (1) it generally requires the use of periodic boundary conditions, in which the simulation describes an infinite threedimensional lattice, with each lattice cell containing a copy of the simulated system, and (2) calculation of the FFT becomes a considerable performance bottleneck to the parallel scalability of MD simulations, due to the manytomany communication pattern employed.
We have developed in NAMD the multilevel summation method (MSM) as an alternative to PME. MSM avoids the use of the FFT in its calculation, instead employing the nested interpolation in real space of softened pair potentials, which permits the use of nonperiodic boundary conditions if desired and enables better parallel scaling with a treelike communication pattern. The implementation of MSM in NAMD is shown to give excellent agreement with PME for many different water properties, especially for interface potentials at airwater and membranewater interfaces for which longrange Coulombic interactions are crucial. Applications demonstrate the suitability of MSM for systems with semiperiodic and nonperiodic boundaries.
Essential Ideas
Figure 1  MSM essential ideas: splitting the Coulomb interaction potential and interpolating the longer range parts from nested grids. 
MSM approximates electrostatic interactions by smoothly splitting the Coulomb interaction potential and interpolating the longrange parts of the potential from grids. The splitting expresses the interaction as a shortrange part plus a series of slowly varying parts of longer range. The shortrange part is calculated exactly between nearby particles. The longrange parts are each interpolated from a grid, where, by construction, the longer ranged parts are interpolated from grids of increased spacing. The interpolation to each atom is a localized operation that requires only nearby grid points. The potential calculated at each grid point requires summing over the charge of the surrounding grid points within the truncation distance for that part of the interaction potential. By nesting the interpolation between grids, the number of terms summed on each level and between levels increases linearly with the number of atoms. The procedure is illustrated in Figure 1. A detailed presentation of the algorithm is available on a separate web page.
Parallel Implementation and Performance
The parallelization of MSM is implemented using Charm++, the parallel language extension of C++ used by NAMD that provides a messagedriven programming paradigm expressed as asynchronously executable objects. NAMD employs a hybrid data and force decomposition for scalable parallelism, decomposing the atoms spatially into uniform patches and also decomposing the work objects for calculating the shortrange nonbonded forces between nearest neighbors. MSM and PME both make use of the existing infrastructure in NAMD for calculating their shortrange parts, modifying the functional form of the interaction to fit each respective method.
Figure 2  Performance of MSM and PME on ApoA1 (92 K atoms) with 1 fs time stepping, on Cray XE6 (32 cores per node). MSM using Hermite interpolation with 5 angstrom grid spacing and 12 angstrom splitting distance achieves 19.6 ns/day with 1536 cores, whereas PME using quartic interpolation with 1.2 angstrom grid spacing and 12 angstrom splitting distance tops out at 17.6 ns/day with 1024 cores. Performance of MSM here exceeds PME at 512 cores and higher. 
The longrange part of MSM is parallelized using a combination of domain decomposition and force decomposition, much like the approach used to parallelize the shortrange part. As with the shortrange force calculation, sufficient work for scalable parallelization is available by carrying out the grid cutoff calculations involving interactions between neighboring blocks of grid points through separately schedulable work objects. The restriction and prolongation calculations introduced above require a much smaller amount of work, so they are best kept on the same processor as their corresponding grid blocks, as the overhead of communicating the data to a remote processor exceeds the runtime of the actual calculation step. Similarly, the anterpolation and interpolation calculations involving atom coordinates are kept processorwise with the patches.
The MSM algorithm is well suited to finegrained parallel execution on CPU SIMD hardware, due to its predictably and uniformly strided memory access patterns, and because its innermost loops over grids can be decomposed into operations on vectors that match the size and stride of the SIMD processing units provided by the underlying computing hardware.
Performance and scaling of MSM and PME are comparable, despite the fact that the NAMD implementation of PME has had many more years of development effort than MSM. The NAMD benchmark system ApoA1 (92 K atoms) simulated with 1 fs single time stepping using varying numbers of Cray XE6 nodes (32 cores per node) of the Blue Waters JYC test computer system shows that, while PME achieves faster simulation rates for smaller core counts, MSM achieves comparable scaling results with slightly faster simulation rates for larger core counts. Figure 2 plots the results showing nanoseconds per day versus number of cores. The performance crossover point occurs between 256 and 512 cores, for which the MSM performance continues to scale up to 1536 cores, whereas the PME performance reaches a plateau at 1024 cores.
Calculating Water Properties and AirWater Interface Potential
Figure 3  Table and plots comparing MSM and PME calculations of various water properties and the airwater interface potential. 
Water is the most fundamental constituent of living cells, making it an excellent test bed for the computational treatment of electrostatic forces. Indeed numerous structural, dynamic, dielectric, and mechanical properties of water can serve to illustrate the accuracy of electrostatic force descriptions by means of MSM and PME, for example, density, radial distribution function, diffusion constant, dielectric constant, distancedependent Kirkwood factor or surface tension. To compare simulation results obtained with MSM and PME electrostatic descriptions, a box of water was simulated for 6 ns in an NPT ensemble with a modified TIP3P water model in the CHARMM force field. The resulting water properties are compared in Figure 3. For all water properties calculated, agreement was found between MSM and PME, indicating that MSM furnishes the same accuracy as PME in electrostatic descriptions.
A key electrostatic property of water arises at the airwater interface in the form of the socalled interface potential. This potential results from an alignment of water dipole and quadrupole moments. Longrange electrostatic interaction is essential for the magnitude of the potential, and only reliable electrostatic treatments can produce it correctly. As shown in the lower righthand plot of Figure 3, MSM gives the same interface potential as PME (in case of either a fully periodic boundary or a semiperiodic boundary), showing that MSM provides a reliable description of the longrange electrostatic interaction.
Applications to Semiperiodic and Nonperiodic Systems
One advantage of MSM over PME is the ability to be applied to systems with semiperiodic or nonperiodic boundary conditions. Introduced here are three exemplary applications for which MSM is more suitable than PME.
Cellular membranes surrounding neurons commonly maintain a voltage gradient for batterylike energy storage. A difference in charge distribution on either side of a membrane, for example, having more cations on one side and anions on the other, creates an electric field, which when strong enough can cause a pore to form in the membrane. Such systems simulated with PME require two different water compartments containing their respective positive and negative ions separated by two membranes, so as to permit a periodic boundary condition along the membrane normal. Rather than doubling the size of the system, MSM can be used without periodicity in the zdirection to simulate a single membrane. Comparing both approaches for a membrane bilayer made of POPC lipids yields the same results: within 1 ns, a water pore develops in the membrane due to the electric field, as shown in Figure 4(a), allowing ions to diffuse through the pore. The MSM approach uses half the number of atoms as the PME approach while describing a more realistic situation.
Figure 4  Applications of MSM to semiperiodic and nonperiodic systems. Image (a) shows a semiperiodic membrane bilayer, for which different ion concentrations above and below the membrane causes a water pore to form and allows ions to diffuse. Image (b) shows a semiperiodic graphene nanopore separating different ion concentrations that diffuse through the nanopore during simulation. Image (c) shows a nonperiodic simulation of a protein in a water droplet. 
Bioengineering is today developing nanoscale sensors for medical diagnostics, for example, for cost effective DNA sequencing. One type of sensor involves nanopores situated in a graphene sheet, the sheet acting much like a cellular membrane. MD simulations offer guidance to nanosensor development by offering engineers microscopic views of the measuring processes involving nanosensors. The ability of MSM to simulate single membranes with semiperiodic boundary conditions can be employed for graphene nanopores, to describe ion conduction through a nanopore embedded into a single graphene sheet. As in the previous example, the sheet separates ions initially at different concentrations above and below the graphene sheet. The image in Figure 4(b) shows the diffusion of ions through the nanopore.
MSM can also be applied to systems without any periodic boundary. For example, proteins can be simulated in a water droplet. The total number of water molecules in a droplet is almost half that required for a cubic box enclosing the sphere, as would be required for simulation with PME. Moreover, simulations with nonzero net charge can be performed with MSM for nonperiodic systems, while the net charge must be chosen to vanish for PMEbased simulations. The image in Figure 4(c) shows a tryptophan cage protein in a water droplet simulated for 1 microsecond using MSM.
Publications

Publications Database Multilevel summation with Bspline interpolation for pairwise interactions in molecular dynamics simulations. David J. Hardy, Matthew A. Wolff, Jianlin Xia, Klaus Schulten, and Robert D. Skeel. Journal of Chemical Physics, 144:114112, 2016. (16 pages). 
Publications Database Multilevel summation method for electrostatic force evaluation. David J. Hardy, Zhe Wu, James C. Phillips, John E. Stone, Robert D. Skeel, and Klaus Schulten. Journal of Chemical Theory and Computation, 11:766779, 2015. 
Publications Database Multilevel summation of electrostatic potentials using graphics processing units. David J. Hardy, John E. Stone, and Klaus Schulten. Journal of Parallel Computing, 35:164177, 2009.
Investigators
Collaborators
Page created and maintained by David Hardy.