hp150111 「京」以外 HPCI 一般利用 HPCI other than K General Use
粒子格子多重極法（PMMM）による大規模分子動力学シミュレーション
Large scale molecular dynamics simulation using Particle Mesh Multipole Method (PMMM)
大野洋介 Yousuke OHNO 理化学研究所 RIKEN 要旨 大規模古典分子動力学シミュレーションにおいて演算量の大半を占める静電相互作用計算を効率 化するアルゴリズムとして粒子格子多重極法（PMMM）[1]が提唱されている。PMMM では高速フ ーリエ変換（FFT）と高速多重極展開法（FMM）を併用することで、FFT の並列化効率と FMM の 演算コストを緩和している。 本研究では PMMM を用いた分子動力学シミュレーションコードを開発し、中規模分子系のシミュ レーションを実施した。現在主流の FFT を使用する粒子格子エバルト（PME）法[2]と比較し、同 等のエネルギー計算精度を得られる条件で PME 法より計算時間が短いことを確認した。 キーワード：粒子格子多重極法、高速フーリエ変換、高速多重極展開、分子動力学、大規模並列 AbstractParticle Mesh Multipole Method (PMMM) [1] is a new algorithm that improves the calculation efficiency of electro-static interaction used in large scale classical molecular dynamics (MD) simulation. The PMMM combines Fast Fourier Transform (FFT) and multipole method and solves the problem of FFT’s parallel efficiency and FMM’s calculation cost. In this work, we developed MD code using PMMM and performed middle-scale MD simulation. Experimental results show that the PMMM requires shorter execution time than Particle Mesh Ewald (PME) method [2] which is a widely used method using FFT, under the condition that the accuracy of the calculated energy is equivalent.
Keywords: Particle Mesh Multipole Method, Fast Fourier Transform, Fast Multipole Method, Molecular
Dynamics, Massive Parallel
1. 研究の背景と目的 ニュートン力学に基づいて原子・分子の多粒子系の時間発展を計算する古典分子動力学 （Molecular Dynamics, 以下 MD）シミュレーションの大規模化では静電相互作用の演算量が問題 となる。静電相互作用は遠距離二体相互作用であるため演算量が粒子数の二乗に比例するからで ある。この演算量の問題を回避するため実際の MD シミュレーションでは遠距離粒子ペアの静電 相互作用演算を削減する近似手法が使用される。静電相互作用は力の大きさが距離の二乗に反比 例する事を利用し、短距離ペアは厳密に直接計算し、遠距離ペアの相互作用は近似計算すること で精度を維持しつつ計算量を削減する。現在主流の手法は高速フーリエ変換（Fast Fourier Transform, 以下 FFT）を用いた Particle Mesh Ewald（PME）法である。PME 法では粒子の電荷を 格子点に割り付け、FFT を使って波数空間に変換することで演算量を削減する。格子点数密度は ほぼ粒子密度に比例するため、FFT の演算量は粒子数 N に対して O(NlogN)となる。演算量は削 減できるが、大規模並列計算においては並列化効率が問題になる。FFT より並列性の高い手法と して高速多重極展開法（FMM）の利用も進みつつある[3]。FMM では粒子の電荷分布を多重極で 表現する。異なる空間スケールで階層的に多重極展開し遠距離相互作用を粗視化することで相互 作用計算量が O(N)となる。相互作用数は削減でき並列環境での通信も FFT より少ないが多重極 展開の計算量が多いため、並列度が小さいと PME 法より計算時間が長くなる。 大中規模分子動力学において並列性と演算量のバランスをとった計算手法として FFT と多重極 法を組み合わせた粒子格子多重極法（PMMM）[1]が提唱されている。PMMM では FMM と同じ ように空間分割領域毎に粒子電荷分布を多重極展開で表現する。この多重極間の相互作用を PME 法のように FFT で波数空間に変換して計算する。この FFT の点数は空間分割領域の個数となる。 また、多重極展開の各次数の項の FFT は独立している。PME 法の場合は空間分割領域内に 1000 点程度の格子点を置くのが典型的であり、PMMM で 5 次の多重極展開を使用する場合項数は 36 個となるため、PME 法と比べると 1/1000 の点数の FFT を 36 個独立に行うことになる。点数が少 ない複数の FFT を独立に行うため PME 法より FFT の並列化が容易となる。階層的な多重極展開 をしないため FMM より多重極の演算が少なくなる。 本研究の目的はPMMMをタンパク質の分子動力学シミュレーションに使用しPME法より高い並 列化効率が得られることを実証することである。
2. 計算モデル
計算対象：水中タンパク質 1337208 原子
先行研究[3]において使用された系で、「京」の強スケール限界を測定できるように粒子数が 選択されている。300K、1 気圧の条件を満たすように体積を変化させる計算を行い、体積は (23.7 nm)3となった。
タンパク質：TTHA1718(PDB ID1 2ROE) 100 原子×8, ovalbumin(PDB ID 1OVA) 5998 原子×64 水(TIP3P: 3 個の粒子で表現した水分子のモデル)×316,952, Na 1200, Cl 480
細胞内環境に近いタンパク質濃度に設定した。 MD 条件
23.7 nm 周期境界、NVE(N:粒子数、V:体積、E:全エネルギーが保存される条件)、時間刻み 1fs/step cutoff (van der Waals force) 1.4 nm
近距離相互作用 計算空間全体を243個に分割した(0.99nm)3の空間分割単位(セル)で管理する。 アルゴリズム個別設定 PME 法 格子点への電荷割付 6 次 B-spline FFT 格子点数 2883 grid 静電相互作用直接計算 cutoff 1.4 nm PMMM パラメータ FMM 次数 5、項数 36、1 セルを 1 個の多重極で展開する（展開空間サイズ 0.99 nm） FFT 243 grid, fftw3 (single process×36)
各パラメータはエネルギーの相対誤差が 10-5程度に収まるように設定した。 古典分子動力学計算コードは、理化学研究所・次世代計算科学研究開発プログラムで開発した cppmd を利用した。近距離相互作用計算、PME 法は cppmd に実装されたものを利用し、PPPM コ ードを追加した。FFT 計算には fftw3 を使用した。 使用した計算機は、平成 27 年度「京」以外の HPCI システム一般課題で提供を受けた、東京大学 情報基盤センターFX10 スーパーコンピュータシステム（Oakleaf-FX）である。プロセッサは SPARC64TM IXfx（動作周波数 1.858 GHz、理論コア性能 14.78 Gflops）、ノードあたり 16 コア、 理論ノード性能 236.5 Gflops、総ノード数 4800 であった。本研究では最大 504 ノードを使用した。
3. 並列計算の方法と効果（性能） 近距離相互作用と遠距離相互作用に独立なプロセスを割り当てた。近距離相互作用プロセスはセ ル単位で 3 次元分割した。 表 1 PME 法プロセス 近距離相互作用 遠距離相互作用 合計 セル数／ノード プロセス数 プロセス数 プロセス数 64 (4×4×4) 216 (6×6×6) 36 252 32 (4×4×2) 432 (6×6×12) 72 504 16 (4×2×2) 864 (6×12×12) 144 1008 表 2 PMMM プロセス 近距離相互作用 遠距離相互作用 合計 セル数／ノード プロセス数 プロセス数 プロセス数 64 (4×4×4) 216 (6×6×6) 36 252 32 (4×4×2) 432 (6×6×12) 36 468 16 (4×2×2) 864 (6×12×12) 36 900 利用した cppmd の仕様により 8 スレッド／プロセスとした。使用したFX10は 1 ノード 16 コア なので 1 ノード 2 プロセスで実行した。
プロセス並列は MPI、スレッド並列は OpenMP で記述した。SIMD 並列はコンパイラによるルー プの SIMD 化機能を使用した。 図 1 に計算時間測定の結果を示した。1 ステップの最短計算時間は、PME 法(8064 コア) 118ms、 PMMM(7200 コア) 70ms で、PMMM は PME 法の 0.59 倍、コア数の違いを考慮すると 0.53 倍で あった。どちらも著しい並列性能の低下は見られず、強スケーリング性能は PME 法 0.79、PMMM 0.73-0.81 であった。 図 1 計算時間：コア数に対する 1 ステップあたりの計算時間（ms） （左：PME 法、右：PMMM）
4. 研究成果 計算精度の検証のため同一の初期条件で長時間計算し、一原子あたりの全エネルギー（運動エネ ルギーとポテンシャルエネルギーの合計）と系全体の平均温度の時間変化をPME 法と比較した。 この MD での温度は原子の運動エネルギーの総和 Ekから T=Ek/(0.5Nk)（N: 系全体の自由度、 k: ボルツマン定数）で計算されるもので、この系では 300K は 1 原子あたりの運動エネルギー 0.640 kcal/mol に相当する。 図 2 エネルギーと温度の時間変化 （左：1 原子あたりの全エネルギー、右：温度の時間変化） 図 2 に最大並列数（PME 法 1008 プロセス、PMMM 900 プロセス）での原子あたりの全エネルギ ーと温度の時間変化を示した。各手法において並列度による違いはなかった。全エネルギーは 10-5 未満の相対誤差で PME 法と一致する一定値であった。温度はどちらも標準偏差 0.2K で変動して いる。1000 ステップ以上の区間平均では PMMM の温度が 0.2K 高い。 平均運動エネルギーに比例する平均温度が一致しないのは個々の粒子の速度が異なるからであ る。ポテンシャルエネルギーにも差があり運動エネルギーの差と相殺されて全エネルギーは一致 している。PPPM と PME 法で力とポテンシャルの誤差が異なるために同一初期条件でも速度と 座標の時間発展に差が生じたと考えられる。 5. まとめと今後の課題 並列化効率に関しては今回の並列度では PME 法の並列化効率低下が見られなかったため PMMM の優位性は確認できなかったが、計算時間は 0.6 倍程度に短縮されており計算時間での優位性は 確認できた。PME 法の並列性能低下が見られる領域での比較が今後の課題である。「京」での PME 法の計算で、1,000,000 原子 3,200 コア(8,000 ノード)で並列効率の低下がみられる実施例もあり[5]、
エネルギー保存性は維持されており絶対値も PME 法と一致しており、今回の PMMM の計算条件 はタンパク質等の分子動力学シミュレーションには十分なものであったと言える。FMM では PME 法と比較してエネルギー保存性が劣っていることが知られており、同一精度を達成するのに 必要な多重極展開の次数や個数は PMMM が有利と考えられるが、定量的な比較は今後の課題で ある。温度は 0.2K の差があり温度揺らぎと同程度であるが定常的な差でありタンパク質の変形 等の動態への影響は今後の課題である。 参考文献
[1] Keigo Nitadori, Particle mesh multipole method: An efficient solver for gravitational/electrostatic forces based on multipole method and fast convolution over a uniform mesh, arXiv:1409.5981v2 [astro-ph.IM]. [2] T. Darden, D. York, L. Pedersen, Particle mesh ewald: An n log(n) method for ewald sums in large
systems, J. Chem. Phys. 98 (1993) 10089-10092.
[3] Yousuke Ohno; Rio Yokota; Hiroshi Koyama; et al., Petascale molecular dynamics simulation using the fast multipole method on K computer, Computer Physics Communications Vol. 185 Issue 10, 2575-2585 (2014).
[4] wwPDB https://www.wwpdb.org/
[5] Chigusa Kobayashi, Jaewoon Jung, Yasuhiro Matsunaga; et al., GENESIS 1.1: A hybrid-parallel molecular dynamics simulator with enhanced sampling algorithms on multiple computational platforms. J. Comput. Chem., 38, 2193-2206 (2017).
