OpenCLで実装されたSPHのパフォーマンス比較

OpenCLでSPH法のプログラムを実装し、様々なプロセッサで性能を比較している。

このSPHのプログラムは、​http://arxiv.org/abs/1112.4539 の論文で発表したGPUでの直接treewalkを拡張したもので、この論文のAppendixにある、OpenCL Cで書かれたコードがベースとなっている。我々の手法は、treewalkだけをGPUなりmulti core CPUで実行するというものであり、我々の意見ではこうするのが一番性能がでる。この論文は重力相互作用場合の性能評価をしているが、この方法の拡張で、SPHの近接相互作用の和を直接計算することができる。粒子をキャッシュにやさしい順番に並び替えるということも同じく効果的である。 さらに、現時点のコードでは、ホスト計算機でのtreeのconstructionを含めて、全体的に実装しなおして高速化した。これは、GPUによりtreewalkが高速になった反面、今度はtree constructionがボトルネックになってきたからである。上記論文のデータでは、80万粒子の場合、tree constructionには約30％の時間がかかっている。アムダールの法則を持ち出すまでもなく、より大きい粒子の計算では、GPUの高速さ(ベクトル計算機として)の意義が失われてしまう。この論文のコードでは、tree constructionには、シングルスレッドで八分木のポインタに基づいた方法を使っていた。これはBarnesのオリジナルのFortranでの実装とほぼ同様である。今のコードではMorton keyを並列ソートしてから、直接Next/Moreのポインタを生成する手法とした(Next,Moreの意味は論文参照)。

このNext/Moreを使う方法は、文献によると"Brother/child representation"と呼ばれるらしい(参照： ​http://onlinelibrary.wiley.com/doi/10.1111/j.1467-8659.2007.01104.x/abstract )。あるいは、CS業界ではfully threaded treeと呼ばれるらしい。treewalkは再帰関数として記述することができ、非常に美しいアルゴリズムの一つであるが、再帰関数はそのまま扱うとスタックの操作が頻繁になるため、ループに比べると性能がよくないかもしれない。末尾再帰はループに変換することができる。つまりNext/Moreのポインタを使うことで、treewalkを再帰を使わずに、また明示的なスタック処理も必要なく扱うことができる。これはGPUに向いたアルゴリズムで(最初に提案されたのはDubinskiの修士論文(1988) or Makino&Hut 1989 ？かもしれない)、Makino(1990)にて提案されたtree法のベクトル化手法の一つである。実はこの時、他に二人の著者(BarnesとHernquist)がそれぞれ別のベクトル化を提案して、この三つの論文はまとめてJournal of Computational Physicsに掲載されている。このMakinoによる末尾再帰のループへの変換の応用といい、そもそもBarnesによる八分木ツリーの提案も、彼らがIASに滞在していた頃におこなわれている（ような気がする）。しかも、これはどちらもConnection Machineに関係していて、そしてそれはMITのLISP/Schemeハッカー達の多大なる影響のもとで発明されたと勝手に思っているが、どうだろうか？機会があれば、本人に尋ねてみよう。リンク。と、ものすごい脱線はさておき。

この方法の利点は、treeのデータ構造として必ずしも八分木を想定しなくてもよいことである。そもそもオリジナルのtree法では、treeの各ノードの粒子数が1個になるまで、再帰的にtreeデータを構築していた。だが、粒子数が増えてきた時には、1個になるまでノードを追加していくのは時間がかかるし、実際問題必要のないことが90年(?)頃に発見されて、特にGRAPEとtree法を組み合わせる場合、通常はノード内の粒子が数十個を下回ったらもうノードを作らないほうが速いし、力の精度もよい。こうすることで、単純には元々のtree法の再帰的動的であるという美しさが失われてしまうが、実際上速いのに越したことはない。我々の方法ではソートされた粒子のリストから、基本的には八分木のデータを生成する。ただし八分木のポインタは使わない。また、下の階層のノードでは、ノード内の粒子数が指定した値を下回った場合には、それ以上の階層を作るのをやめて、代わりにNextポインタで粒子を直接リンクする。既にソートされているものをポインタを使ってリンクしなおすことで、実は効率が落ちてしまうが、こうするとtreewalkカーネルは修正をする必要がない。これは、実効的には、最下層ノードはN分木になっていて、その上は八分木になっているようなツリーになる。またポインタを使うことで、上の階層の八分木も保持する必然はなくなる。幾何学的に局所性が保たれるのであればN分木に変更してもよいということになる。これが意味を持つのは並列化をする時である。

さて、性能評価は以下の図。図は随時アップデートされるかもしれない。やっている計算は、いわゆる冷たいガス球の収縮で、hydrodyamicsとself gravityを解いている。状態方程式は理想気体で、もちろん断熱である。粒子数を横軸として1 stepあたりの平均計算時間をプロットした。計算時間はメインループの計算時間より算出しているので、当然tree constructionやGPUとのデータ転送等すべてを含む。演算は全て単精度。IntelのCPUではIntelのSDKを、AMDのGPUとCPUはAMDのSDKを、NVIDIAのGPUではNVIDIAのSDKを使っている。

以下コメント：

• OpenCLのSDKは使いものになる。この結果は、コードベースはひとつで各社のSDKで問題なく動作した。
• GPUではTahitiが一番速い。
• FermiのほうがGT200より速い。
• 24 core機がCPUの中では一番速い。
• Sandybridge-Eは24 coreのOpteronとあまり性能が変わらない。
• APUの結果はGPUを使った場合。Nの小さい所で大健闘している。次世代のAPUは面白いかも。

。。。

コメントはこちら ​https://plus.google.com/117333869988556483969/posts/Y144K5p53LM