^{1,a)}, Tetsu Narumi

^{2}and Kenji Yasuoka

^{1}

### Abstract

An IPS/Tree method which is a combination of the isotropic periodic sum (IPS) method and tree-based method was developed for large-scale molecular dynamics simulations, such as biological and polymer systems, that need hundreds of thousands of molecules. The tree-based method uses a hierarchical tree structure to reduce the calculation cost of long-range interactions. IPS/Tree is an efficient method like IPS/DFFT, which is a combination of the IPS method and FFT in calculating large-scale systems that require massively parallel computers. The IPS method has two different versions: IPSn and IPSp. The basic idea is the same expect for the fact that the IPSn method is applied to calculations for point charges, while the IPSp method is used to calculate polar molecules. The concept of the IPS/Tree method is available for both IPSn and IPSp as IPSn/Tree and IPSp/Tree. Even though the accuracy of the Coulomb forces with tree-based method is well known, the accuracy for the combination of the IPS and tree-based methods is unclear. Therefore, in order to evaluate the accuracy of the IPS/Tree method, we performed molecular dynamics simulations for 32 000 bulk water molecules, which contains around 10^{5} point charges. IPSn/Tree and IPSp/Tree were both applied to study the interaction calculations of Coulombic forces. The accuracy of the Coulombic forces and other physical properties of bulk water systems were evaluated. The IPSp/Tree method not only has reasonably small error in estimating Coulombic forces but the error was almost the same as the theoretical error of the ordinary tree-based method. These facts show that the algorithm of the tree-based method can be successfully applied to the IPSp method. On the other hand, the IPSn/Tree has a relatively large error, which seems to have been derived from the interaction treatment of the original IPSn method. The self-diffusion and radial distribution functions of water were calculated each by both the IPSn/Tree and IPSp/Tree methods, where both methods showed reasonable agreement with the Ewald method. In conclusion, the IPSp/Tree method is a potentially fast and sufficiently accurate technique for predicting transport coefficients and liquid structures of water in a homogeneous system.

K.T. was supported by Grant-in-Aid for Japan Society for the Promotion of Science (JSPS) Fellows 21-5452 of the Ministry of Education, Culture, Sports, Science, and Technology (MEXT). T.N. and K.Y. were supported by the Core Research for the Evolution Science and Technology (CREST) of the Japan Science and Technology Corporation (JST).

I. INTRODUCTION

II. METHODS

A. Tree-based method

B. IPS method

C. IPS/Tree method

D. MD simulations

III. RESULTS AND DISCUSSION

A. Accuracy of Coulombic forces

B. Potential energy and self-diffusion

C. Radial distribution function

IV. CONCLUSION

### Key Topics

- Molecular dynamics
- 14.0
- Self diffusion
- 10.0
- Cell division
- 5.0
- Computer simulation
- 5.0
- Particle distribution functions
- 5.0

## Figures

An example of the hierarchical cell for the tree-based method. The root level cell (the simulation cell) is recursively sub-divided into octants. Each *P* cell is a level 2 cell, which is divided twice from the root level cell. Cell *A*, shaded cells, and *M* cells are level 3 cells. In this figure, level 3 cells mean leaf level cells. Shaded cells are nearest neighbor cells, and *M* cells are the next nearest neighbor cells from cell *A*.

An example of the hierarchical cell for the tree-based method. The root level cell (the simulation cell) is recursively sub-divided into octants. Each *P* cell is a level 2 cell, which is divided twice from the root level cell. Cell *A*, shaded cells, and *M* cells are level 3 cells. In this figure, level 3 cells mean leaf level cells. Shaded cells are nearest neighbor cells, and *M* cells are the next nearest neighbor cells from cell *A*.

(a) In the IPS/Tree method, the simulation cell is recursively sub-divided, like the typical tree-based method. Shaded cells are local region cells that partially include the cutoff radius territory. Cells *M* and *P* are both middle-range cells. The scheme of the interaction calculation between *A* cells and middle-range cells is always the same regardless of the difference of the size of each middle-range cells. (b) All positive point charges inside each divided cell was plotted. If the multipole expansion is not applied, interaction calculations between cells *A* and numerous point charges are needed. (c) All positive point charges after the P^{2}M^{2} application are plotted. The multipole expansion up to the pth order can theoretically be reproduced by the expansion of *K* _{min}(*p*) = ⌈(*p* + 1)^{2}/4⌉ pseudoparticles. For *p* = 2, the minimum number of pseudoparticles necessary is *K* _{min}(2) = 3 for each cell. (d) The IPS method is applied for point charges after the P^{2}M^{2} application, which includes pseudoparticles. The corner regions (outside of the large cutoff radius) is an isotropic periodic force field. Since the pseudoparticle is a mere point charge, interaction calculations using approximated IPS potentials by polynomial functions can be easily performed.

(a) In the IPS/Tree method, the simulation cell is recursively sub-divided, like the typical tree-based method. Shaded cells are local region cells that partially include the cutoff radius territory. Cells *M* and *P* are both middle-range cells. The scheme of the interaction calculation between *A* cells and middle-range cells is always the same regardless of the difference of the size of each middle-range cells. (b) All positive point charges inside each divided cell was plotted. If the multipole expansion is not applied, interaction calculations between cells *A* and numerous point charges are needed. (c) All positive point charges after the P^{2}M^{2} application are plotted. The multipole expansion up to the pth order can theoretically be reproduced by the expansion of *K* _{min}(*p*) = ⌈(*p* + 1)^{2}/4⌉ pseudoparticles. For *p* = 2, the minimum number of pseudoparticles necessary is *K* _{min}(2) = 3 for each cell. (d) The IPS method is applied for point charges after the P^{2}M^{2} application, which includes pseudoparticles. The corner regions (outside of the large cutoff radius) is an isotropic periodic force field. Since the pseudoparticle is a mere point charge, interaction calculations using approximated IPS potentials by polynomial functions can be easily performed.

(a) θ dependence of *e* _{Ewald}. IPSp/Tree at *r* _{c} ⩾ 2.7 nm has nearly the same error as IPSp with a long cutoff radii in comparison with the Ewald sum. On the other hand, the IPSn/Tree gives a relatively large error for estimating the Coulombic interaction. This is from the difference of the interaction treatment between IPSn and IPSp. (b) θ dependence of *e* _{IPS}. The interaction of the typical tree-based method with second-order multipole expansion has 10^{−2} − 10^{−4} order relative error, and this error is roughly proportional to θ^{−3}. The error of IPSp/Tree has almost the same dependence against θ. In the case of IPSn/Tree, however, 10^{−3} order error remains even if the short cutoff radius has a value larger than 2.7 nm, so θ dependence of the error is weaker than that of the IPSp/Tree and other typical tree-based methods.

(a) θ dependence of *e* _{Ewald}. IPSp/Tree at *r* _{c} ⩾ 2.7 nm has nearly the same error as IPSp with a long cutoff radii in comparison with the Ewald sum. On the other hand, the IPSn/Tree gives a relatively large error for estimating the Coulombic interaction. This is from the difference of the interaction treatment between IPSn and IPSp. (b) θ dependence of *e* _{IPS}. The interaction of the typical tree-based method with second-order multipole expansion has 10^{−2} − 10^{−4} order relative error, and this error is roughly proportional to θ^{−3}. The error of IPSp/Tree has almost the same dependence against θ. In the case of IPSn/Tree, however, 10^{−3} order error remains even if the short cutoff radius has a value larger than 2.7 nm, so θ dependence of the error is weaker than that of the IPSp/Tree and other typical tree-based methods.

The results of oxygen-oxygen, oxygen-hydrogen, and hydrogen-hydrogen radial distribution functions for (a) IPSn/Tree and (b) IPSp/Tree method at *r* _{c} = 1.5 nm. All radial distribution functions sufficiently follow the results of the Ewald sum.

The results of oxygen-oxygen, oxygen-hydrogen, and hydrogen-hydrogen radial distribution functions for (a) IPSn/Tree and (b) IPSp/Tree method at *r* _{c} = 1.5 nm. All radial distribution functions sufficiently follow the results of the Ewald sum.

The magnified figures of oxygen-oxygen, oxygen-hydrogen, and hydrogen-hydrogen radial distribution functions for (a) IPSn/Tree and (b) IPSp/Tree method at *r* _{c} = 1.5 nm. The detail of the shape of RDFs at 0.86 nm ⩽ *r* ⩽ 5.0 nm was plotted. There were no special deviation of *g* _{IPSn/Tree}(*r*) and *g* _{IPSp/Tree}(*r*) around *r* _{c} and *R* _{c}. The IPSn/Tree can avoid innate IPSn problems where a sharp deviation appears.^{25,27,28}

The magnified figures of oxygen-oxygen, oxygen-hydrogen, and hydrogen-hydrogen radial distribution functions for (a) IPSn/Tree and (b) IPSp/Tree method at *r* _{c} = 1.5 nm. The detail of the shape of RDFs at 0.86 nm ⩽ *r* ⩽ 5.0 nm was plotted. There were no special deviation of *g* _{IPSn/Tree}(*r*) and *g* _{IPSp/Tree}(*r*) around *r* _{c} and *R* _{c}. The IPSn/Tree can avoid innate IPSn problems where a sharp deviation appears.^{25,27,28}

## Tables

Parameters of the Coulombic potential for the IPSn method. The average deviation of the polynomial fitting function, Eq. (7), from the analytic solution Eq. (8), is about 2 × 10^{−8} *q* _{ i } *q* _{ j }/*R* _{c}.

Parameters of the Coulombic potential for the IPSn method. The average deviation of the polynomial fitting function, Eq. (7), from the analytic solution Eq. (8), is about 2 × 10^{−8} *q* _{ i } *q* _{ j }/*R* _{c}.

Potential energy and self-diffusion coefficient for the IPS/Tree method, IPS method with long cutoff, and Ewald sum.

Potential energy and self-diffusion coefficient for the IPS/Tree method, IPS method with long cutoff, and Ewald sum.

Article metrics loading...

Full text loading...

Commenting has been disabled for this content