GPUでsingle spin flip

single spin flipによる状態更新はスピン系のモンテカルロ計算で基本中の基本ですが,
GPUでsingle spin flipによる状態更新をした場合の情報が見つからなかったので,
自分でプログラムを書いて性能の比較をしてみました.

性能を比較するにあたり, どのスピン模型を使うか迷ったのですが, 2次元の強磁性XY模型
で比較することにしました. 「標準」という観点から模型を選ぶとIsing模型かな, とも考
えたのですが, Ising模型の場合, bit codingなどの高速化も可能なわけで, 何をもって正当な
比較と言えるのかよくわからなくなったので, そういった側面も考慮しての模型の選択と
なりました.

比較対象となるCPUでの計算ですが, 正方格子をAとBの2つの副格子に分け, OpenMPの
8並列で実行速度を計算しました. 乱数はdSFMT ver2.2.1を使っていますが, elapsed timeに
占める割合はほぼゼロと考えていただいて差し支えありません.

　CPU: intel Core i7 3820 (4core)
　システムサイズ: 1024×1024
　ステップ数: 10000
　random number generator: dSFMT ver2.2.1 (MT19937)
　elapsed time: 246.2 s (= 23.5 ns/spin)

GPUでの計算では, i番目の状態更新時に参照したスピン変数のうち, i+1番目でも利用
するものをshared memoryに格納する工夫だけ行いました. 乱数はMTGP ver1.1.2を
使っていますが, こちらもCPU同様, elapsed timeに占める割合はほぼゼロになっています.

　GPU: NVIDIA GTX680
　システムサイズ: 1024×1024
　ステップ数: 10000
　random number generator: MTGP ver1.1.2 (MT23209)

elapsed timeは1ブロックあたりのスレッド数(行)と総ブロック数(列)の組み合わせに依存します.

	512	256	128	64	32	16
32	15.3	15.8	16.2	21.1	36.0	---
64	---	13.9	14.6	16.2	21.0	35.3
128	---	---	14.0	15.0	16.1	22.1
256	---	---	---	14.6	15.2	17.7
512	---	---	---	---	16.9	16.9

(表中の数字の単位は秒)

なので, この表中でもっとも良い組み合わせはブロック当たりのスレッド数が64で, 総ブロック数が
256ということになり, elapsed timeは13.9秒ということになります. しかし, 13.9秒の中に
は, 実行までのオーバーヘッドがそれなりの割合で含まれていて, スレッド数64, 総ブロック数256
の組み合わせで100000ステップまで計算すると, elapsed timeは124.9秒で10倍より少し小さく
なります. まとめると,

　GPU: NVIDIA GTX680
　システムサイズ: 1024×1024
　ステップ数: 100000
　実行条件: 256 Blocks, 64 threads/Block
　random number generator: MTGP ver1.1.2 (MT23209)
　elapsed time: 124.9 s (= 1.19 ns/spin)

1スピンあたりの更新時間でCPUとGPUの結果の比をとると,

$\displaystyle \frac{23.5}{1.19} \approx 19.7$

で, GPUを利用することで約20倍の高速化を達成することができました.
CPUやGPUで実行したプログラムに手を入れて, 改良することで, 上の
数字は変わり得ますが, CPUとGPUでどちらもあまりがんばっていない
プログラムでの結果として参考にしていただければと思います.