Cell Multipole Method Research Articles

A recursive cell multipole method for atomistic electrodynamics models.

For large plasmonic nanoparticles, retardation effects become important once their length becomes comparable to the wavelength of light. However, most models do not incorporate retardation effects due to the high computational cost of solving for the optical properties of large atomistic electrodynamics systems. In this work, we derive and implement a recursive fast multipole method (FMM) in Cartesian coordinates that includes retardation effects. In this method, higher-order electrodynamic interaction tensors used for the FMM are calculated recursively, thus greatly reducing the implementation complexity of the model. This method allows for solving of the optical properties of large atomistic nanoparticles with controlled accuracy; in practice, taking the expansion to the fifth order provides a good balance of accuracy and computational time. Finally, we study the effects retardation has on the near- and far-field properties of large plasmonic nanoparticles with over a million atoms using this method. We specifically focus on nanorods and their dimers, which are known to generate highly confined fields in their junctions. In the future, this method can be applied to simulations in which accurate near-field properties are required, such as surface-enhanced Raman scattering.

Molecular dynamics study of pattern transfer in nanoimprint lithography

A molecular dynamics simulations model of nanoimprint lithography (NIL) is proposed in order to study the pattern transfer and its related phenomena. The proposed model is similar to a real NIL process imprinting an α-quartz stamp with a rectangular line pattern into an amorphous poly-(methylmethacrylate) (PMMA) film. The polymer deformation behavior and the adhesion and friction effects between the stamp and the polymer film are investigated and their dependency on the pattern aspect ratio is discussed. Force fields including bond, angle, torsion, van der Waals, and electrostatic potentials are used to describe intermolecular and intramolecular interacting forces. Nose-Hoover thermostat is used to control the temperature of the polymer film and cell multipole method is adopted to treat long range interactions. The deformation of the polymer film is observed for two stamps having different aspect ratio patterns. The distributions of density and stress in the polymer film are calculated for the detail analysis of deformation behavior. For a high aspect ratio pattern (aspect ratio = 2.5, imprint depth = 8.0 nm), large amount of springback of the residual polymer film is observed, which is mainly due to the residual compressive stress left in the polymer film. However, for a low aspect ratio pattern (aspect ratio = 1.0, imprint depth = 3.0 nm), the springback is not observed. In addition, adhesion and friction forces are obtained by dividing the polymer film into subregions and calculating the interacting force between each subregion and the stamp. While the adhesion force is nearly constant regardless of the pattern aspect ratio, the friction force increases as the pattern aspect ratio grows, so the friction force becomes larger than the adhesion force when the pattern aspect ratio increases.

Cell multipole method for molecular simulations in bulk and confined systems

One of the bottlenecks in molecular simulations is to treat large systems involving electrostatic interactions. Computational time in conventional molecular simulation methods scales with O(N2), where N is the number of atoms. With the emergence of new simulations methodologies, such as the cell multipole method (CMM), and massively parallel supercomputers, simulations of 10-million atoms or more have been performed. In this work, the optimal hierarchical cell level and the algorithm for Taylor expansion were recommended for fast and efficient molecular dynamics simulations of three-dimensional (3D) systems. CMM was then extended to treat quasi-two-dimensional (2D) systems, which is very important for condensed matter physics problems. In addition, CMM was applied to grand canonical ensemble Monte Carlo simulations for both 3D and 2D systems. Under the optimal conditions, our results show that computational time is approximately linear with N for large systems, average error in total potential energy is about 0.05% for 3D and 0.32% for 2D systems, and the RMS force error is 0.27% for 3D and 0.43% for 2D systems when compared with the Ewald summation.

A large scale molecular dynamics simulation code using the fast multipole algorithm (FMD): performance and application

We present the performance of the fast classical molecular dynamics (MD) code, fast molecular dynamics (FMD), designed for efficient, object-oriented, and scalable large scale simulations, and summarize its application to a liquid crystalline cluster. FMD uses an implementation of the three-dimensional fast multipole method, developed in our group. The fast multipole method offers an efficient way (order O( N)) to handle long range electrostatic interactions, thus, enabling more realistic simulations of large molecular systems. Performance testing was carried out on IBM SP2, SGI Origin 2000, and CRAY T3E massively parallel systems using the MPI massage passing library. The electrostatic forces were tested on models of up to 100,000 randomly placed charges, and on protein and liquid crystalline molecular systems of over 99,000 atoms. Tests on the stability of the method are presented, along with comparisons with direct calculations, the namd2 code, and the physical multipole-based cell-multipole method.

A cell multipole based domain decomposition algorithm for molecular dynamics simulation of systems of arbitrary shape

A domain decomposition algorithm for molecular dynamics simulation of atomic and molecular systems with arbitrary shape and non-periodic boundary conditions is described. The molecular dynamics program uses cell multipole method for efficient calculation of long range electrostatic interactions and a multiple time step method to facilitate bigger time steps. The system is enclosed in a cube and the cube is divided into a hierarchy of cells. The deepest level cells are assigned to processors such that each processor has contiguous cells and static load balancing is achieved by redistributing the cells so that each processor has approximately same number of atoms. The resulting domains have irregular shape and may have more than 26 neighbors. Atoms constituting bond angles and torsion angles may straddle more than two processors. An efficient strategy is devised for initial assignment and subsequent reassignment of such multiple-atom potentials to processors. At each step, computation is overlapped with communication greatly reducing the effect of communication overhead on parallel performance. The algorithm is tested on a spherical cluster of water molecules, a hexasaccharide and an enzyme both solvated by a spherical cluster of water molecules. In each case a spherical boundary containing oxygen atoms with only repulsive interactions is used to prevent evaporation of water molecules. The algorithm shows excellent parallel efficiency even for small number of cells/atoms per processor.

Time-Trimming Tricks for Dynamic Simulations: Splitting Force Updates to Reduce Computational Work

Time-Trimming Tricks for Dynamic Simulations: Splitting Force Updates to Reduce Computational Work

Efficient parallel algorithms for molecular dynamics simulations using variable charge transfer electrostatic potentials

We present an algorithm developed for the efficient computation on parallel computers of molecular dynamics simulations of materials with electrostatic potentials that allow for dynamic charge evolution. This algorithm combines a hierarchical cell multipole method for the fast evaluation of the electrostatic potential with a Broyden–Fletcher–Goldfarb–Shanno scheme for approximate solutions to arbitrary precision of the electronegativity equalization condition, which defines charge transfer in the material. We apply this algorithm to the simulation of a model α-alumina slab at equilibrium. First, we demonstrate that simple fixed-charge models yield different bulk energies than do dynamic-charge models. Second, we compare our algorithm with two previously published techniques for the simulation of dynamic charges in electrostatic materials, discussing both computational speed and numerical accuracy. We show the relative computational speed-up of this scheme over a direct evaluation of the linear equations, and show that intrinsic flaws in other approximate dynamic-charge schemes can be avoided with this approach. © 2000 John Wiley & Sons, Inc. Int J Quant Chem 80: 733–742, 2000

Efficient parallel algorithms for molecular dynamics simulations using variable charge transfer electrostatic potentials

We present an algorithm developed for the efficient computation on parallel computers of molecular dynamics simulations of materials with electrostatic potentials that allow for dynamic charge evolution. This algorithm combines a hierarchical cell multipole method for the fast evaluation of the electrostatic potential with a Broyden–Fletcher–Goldfarb–Shanno scheme for approximate solutions to arbitrary precision of the electronegativity equalization condition, which defines charge transfer in the material. We apply this algorithm to the simulation of a model α-alumina slab at equilibrium. First, we demonstrate that simple fixed-charge models yield different bulk energies than do dynamic-charge models. Second, we compare our algorithm with two previously published techniques for the simulation of dynamic charges in electrostatic materials, discussing both computational speed and numerical accuracy. We show the relative computational speed-up of this scheme over a direct evaluation of the linear equations, and show that intrinsic flaws in other approximate dynamic-charge schemes can be avoided with this approach. © 2000 John Wiley & Sons, Inc. Int J Quant Chem 80: 733–742, 2000

Simulation of a PEO based solid polyelectrolyte, comparison of the CMM and the Ewald summation method

Modeling of polyelectrolytes is very time consuming because of the strong long-range interactions between ions. Traditionally the Ewald summation method is used to calculate the non-bonding Van der Waals and Coulombic interactions. Recently a new much faster method called the cell multipole method (CMM), has been developed to model the non-bonding interactions. The use of the CMM makes the calculations dramatically faster compared to the duration of the calculation when the Ewald summation method is used. In this paper we compare the results of a spectral analysis, pair correlation study and dihedral angle study made for a solid amorphous polyelectrolyte by using both the Ewald summation method and the CMM. The simulated system contained poly(ethylene oxide) sulfonic acid anions, poly(ethylene oxide) (PEO), cations and water. There was no remarkable difference between the results obtained by using Ewald summation method or the CMM. Thus in future the calculations of similar systems can be made applying the fast CMM.

Localized-Density-Matrix Method and Its Application to Carbon Nanotubes

The localized-density-matrix method (Yokojima, S.; Chen, G. H. Chem. Phys. Lett. 1998, 292, 379) is employed to simulate the optical responses of very large carbon nanotubes and polyacetylene oligomers containing 10 000 carbon atoms. The Pariser−Parr−Pople Hamiltonian is used to describe the π electrons in these systems, and the time-dependent Hartree−Fock approximation is employed to calculate the linear optical responses. In the calculation, the fast multipole method or the cell multipole method is employed to evaluate the effects of Coulomb interaction. It is illustrated that the computational time scales linearly with the system size for carbon nanotubes while high accuracy is achieved.

Improved protein free energy calculation by more accurate treatment of nonbonded energy: Application to chymotrypsin inhibitor 2, V57A

We developed a software package for improved free energy calculation, in which spherical solvent boundary potential, cell multipole method, and Nosé-Hoover equation are employed. The performance of the developed software package is demonstrated in the case of valine to alanine mutation of the 57th residue in chymotrypsin inhibitor 2. By using this package, we obtained results quantitatively comparable to experimental results. By the free energy component analysis, it is shown that leucine 51, arginine 65, arginine 67, and phenylalanine 69 residues contribute significantly to the total free energy shift, ΔΔG. Among them, contribution from the hydrophilic arginine 67 residue, which is in close contact with the mutation site, is the largest. Structure around the mutation site is largely changed by the mutation. The structure change is caused mainly by two effects, hydrophobic interaction and short-range interaction along the sequence. Effects of Nosé-Hoover algorithm and Kirkwood reaction field are also discussed. Proteins 30:388–400, 1998. © 1998 Wiley-Liss, Inc.

Improved protein free energy calculation by more accurate treatment of nonbonded energy: Application to chymotrypsin inhibitor 2, V57A

We developed a software package for improved free energy calculation, in which spherical solvent boundary potential, cell multipole method, and Nosé-Hoover equation are employed. The performance of the developed software package is demonstrated in the case of valine to alanine mutation of the 57th residue in chymotrypsin inhibitor 2. By using this package, we obtained results quantitatively comparable to experimental results. By the free energy component analysis, it is shown that leucine 51, arginine 65, arginine 67, and phenylalanine 69 residues contribute significantly to the total free energy shift, deltadeltaG. Among them, contribution from the hydrophilic arginine 67 residue, which is in close contact with the mutation site, is the largest. Structure around the mutation site is largely changed by the mutation. The structure change is caused mainly by two effects, hydrophobic interaction and short-range interaction along the sequence. Effects of Nosé-Hoover algorithm and Kirkwood reaction field are also discussed.

Molecular dynamics for very large systems on massively parallel computers: The MPSim program

We describe the implementation of the cell multipole method CMM in a . . complete molecular dynamics MD simulation program MPSim for massively parallel supercomputers. Tests are made of how the program scales with size . . linearly and with number of CPUs nearly linearly in applications involving up to 10 7 particles and up to 500 CPUs. Applications include estimating the surface tension of Ar and calculating the structure of rhinovirus 14 without requiring icosahedral symmetry. Q 1997 by John Wiley & Sons, Inc.

Implementing the cell multipole method for dipolar and charged dipolar systems

Recently, substantial progress has been made in reducing the cost of computing Coulomb sums in fields such as astrophysics and molecular dynamics (MD). This has culminated in the fast multipole method (FMM), and its Cartesian counterpart, the cell multipole method (CMM). In MD simulations, these methods provide efficient and accurate algorithms for computing the energy and forces of charge-charge interactions. However, dipolar interactions are as important and as challenging to compute. In this paper we implement the cell multipole method (CMM) for dipolar domains. This CMM implementation is an accurate and efficient approach for computing interactions in million-particle dipolar systems. In addition, we use CMM to treat systems of permanent dipoles and induced dipoles, permanent charges and permanent dipoles, and permanent charges and induced dipoles. We derive and discuss expressions for accurate and efficient computation of all interactions in these systems. The approach is also applicable to polarizable charged dipolar systems. Finally, we apply one-level CMM and hierarchical CMM to arbitrary dipolar systems and discuss accuracy and speed. The hierarchical implementation of CMM is shown to scale linearly with the number of dipoles.

Efficient dipole iteration in polarizable charged systems using the cell multipole method and application to polarizable water

The fast multipole method (FMM) and its Cartesian version, the cell multipole method (CMM), allow efficient computation of Coulomb energy and forces in molecular dynamics (MD) simulations. In another work [R. Kutteh and J.B. Nicholas, Implementing the cell multipole method for dipolar and charged dipolar systems, Comput. Phys. Commun., this issue], we implemented CMM for efficient and accurate computation of the energy and forces in dipolar and charged dipolar systems. While CMM provides an efficient way of computing energy and forces in large polarizable systems, the self-consistent dipole iteration remains a computational bottleneck in simulations of these systems. In this paper, we apply CMM (single-level) to substantially reduce the single-iteration time. We further reduce the single-iteration time using a refined CMM iterative algorithm. In addition, we introduce an iterative predictor scheme to reduce the number of iterations to self-consistency. We apply these algorithms to polarizable water using the polarizable SPC (PSPC) model. The “refined CMM/predictor” iterative algorithm is shown to be accurate and extremely efficient. It removes the dipole iteration computational bottleneck in simulations of large polarizable systems.

A generalized fast multipole approach for Hartree—Fock and density functional computations

The fast multipole method (FMM) and the cell multipole method (CMM) compute charge—charge interactions in a system of N particles with O ( N) scaling. We propose an approach for incorporating CMM or FMM into Hartree—Fock (HF) and density functional theory (DFT) calculations. We show how the rate limiting computation of the Coulomb term in HF and DFT can be efficiently carried out using a generalized fast multipole algorithm applicable to higher-order point multipoles.

Multiresolution molecular dynamics algorithm for realistic materials modeling on parallel computers

For realistic modeling of materials, a molecular-dynamics (MD) algorithm is developed based on multiresolutions in both space and time. Materials of interest are characterized by the long-range Coulomb, steric and charge-dipole interactions as well as three-body covalent potentials. The long-range Coulomb interaction is computed with the fast multipole method. For bulk systems with periodic boundary conditions, infinite summation over repeated image charges is carried out with the reduced cell multipole method. Short- and medium-range non-Coulombic interactions are computed with the multiple time-step approach. A separable tensor decomposition scheme is used to compute three-body potentials. For a 4.2 million-particle SiO 2 system, one MD step takes only 4.8 seconds on the 512-node Intel Touchstone Delta machine and 10.3 seconds on 64 nodes of an IBM SP1 system. The constant-grain parallel efficiency of the program is η′ = 0.92 and the communication overhead is 8% on the Delta machine. On the SP1 system, η′ = 0.91 and communication overhead is 7%.

Protein simulations using techniques suitable for very large systems: The cell multipole method for nonbond interactions and the Newton‐Euler inverse mass operator method for internal coordinate dynamics

Two new methods developed for molecular dynamics simulations of very large proteins are applied to a series of proteins ranging up to the protein capsid of tomato bushy stunt virus (TBSV). For molecular dynamics of very large proteins and polymers, it is useful to carry out the dynamics using internal coordinates (say, torsions only) rather than Cartesian coordinates. This allows larger time steps, eliminates problems with the classical description of high energy modes, and focuses on the important degrees of freedom. The resulting equation of motion has the form. [formula: see text] where for T is the vector of generalized forces, M(theta) is the moments of inertia tensor, theta is the vector of torsions, and C is a vector containing Coriolis forces and nonbond forces. The problem is that to calculate the acceleration vector theta from M, C, and T requires inverting M(theta), an order N3 calculation. Since the number of degrees of freedom might be 300,000 for a million atom system, solving these equations every time step is impractical, restricting internal coordinate methods to small systems. The new method, Newton-Euler Inverse Mass Operator (NEIMO) dynamics, constructs the torsional accelerations vector theta = M-1 (T-C) directly by an order N process, allowing internal-coordinate dynamics to be solved for super larger (million atom) systems. The first use of the NEIMO method for molecular dynamics of proteins is presented here. A second serious difficulty for large proteins is calculation of the nonbond forces. We report here the first application to proteins of the new Cell Multipole Method (CMM) to evaluate the Coulomb and van der Waals interactions. The costs of CMM scales linearly with the number of particles while retaining an accuracy significantly better than standard nonbond methods (involving cutoffs). Results for NEIMO and CMM are given for simulations of a wide range of peptide and protein systems, including the protein capsid of TBSV with 488,000 atoms. The computational times for NEIMO and CMM are demonstrated to scale linearly with size. With NEIMO the dynamics time steps can be as large as 20 fs (for small peptides), much larger than possible with standard Cartesian coordinate dynamics. For TBSV we considered both the normal form and the high pH form, in which the Ca2+ ions are removed. These calculations lead to a contraction of the protein for both forms (probably because of ignoring the RNA core not observed in the X-ray).

Atomic level simulations on a million particles: The cell multipole method for Coulomb and London nonbond interactions

The N2 computations implicit in the Coulomb and other long range interactions remain the critical bottleneck in atomic-level simulations of the structure and dynamics of large systems. We report here the cell multipole method which scales linearly with N and requires only modest memory. To demonstrate the feasibility of this approach, we report systematic calculations on realistic polymer systems with up to 1.2 million atoms on a laboratory workstation. The method becomes faster than the exact method for systems of 300 atoms, and for a 1.2 million-atom polymer, it is 2377 times faster. The method treats a class of interactions of the form qiqj/rijp, which includes Coulomb (p=1), London dispersion (p=6), or shielded Coulomb (p=2) interactions. This method is well suited for highly parallel and vector computers.

The reduced cell multipole method for Coulomb interactions in periodic systems with million-atom unit cells

Standard methods for calculating Coulomb interactions of periodic systems use Ewald-type formulations or minimum image approximations, neither of which is practical for large (million-atom) systems. We describe the reduced cell multipole method which is 38 times faster than the Ewald method for a 4816-atom (unit cell) system. The method scales linearly with the number of atoms in the cell, making simulations of million-atom systems practical, as demonstrated in a series of applications on a laboratory workstation.