はじめに
そうとうがんばってSIMD化して、もう速くならないだろうと思ってたらkohnakagawaさんから「もう少し高速化しました」というプルリクが来て慌てた話。
コードは
https://github.com/kaityo256/lj_simdstep
においてある。
どこを直したか
一度に256ビット取ってくる命令を使うため、配列をパディングして、(x,y,z,0)という並びにしていた。それを4要素のベクトルとしてロードして、相対座標ベクトルを作るのだが、ループを4倍展開しているので、
(dx_1, dy_1, dz_1, 0)\\
(dx_2, dy_2, dz_2, 0)\\
(dx_3, dy_3, dz_3, 0)\\
(dx_4, dy_4, dz_4, 0)
という四本のベクトルができる。これを、それぞれ自乗してからシャッフルして足すことで、相対距離を作っていた。
(dr^2_1, dr^2_1, dr^2_1,0) \\
(dr^2_2, dr^2_2, dr^2_2,0) \\
(dr^2_3, dr^2_3, dr^2_3,0) \\
(dr^2_4, dr^2_4, dr^2_4,0) \\
で、これを転置して、
v_r = (dr^2_1, dr^2_2, dr^2_3, dr^4)
を作ってから力の計算を4並列でやっていた。しかし、どうせ最後に転置するなら、最初に転置した方が計算量が減る。つまり、
(dx_1, dy_1, dz_1, 0)\\
(dx_2, dy_2, dz_2, 0)\\
(dx_3, dy_3, dz_3, 0)\\
(dx_4, dy_4, dz_4, 0)
を転置して、
v_x = (dx_1, dx_2, dx_3, dx_4)\\
v_y = (dy_1, dy_2, dy_3, dy_4)\\
v_z = (dz_1, dz_2, dz_3, dy_4)\\
v_0 = (0, 0, 0, 0)
を得たら、あとは
v_r = v_x^2+v_y^2+v_z^2
でできるじゃん、というもの。言われてみればそれはそうだ。なんで気が付かなかったんだろ○|￣|＿
結果
というわけで上記を実装したものをこれまでと同条件で測定。
| 計算手法 | 実行時間[sec] | 備考 |
|---|---|---|
| pair | 8.201647 | ペアごとにループをまわしたもの |
| pair_swp | 7.094673 | 上記をSWPしたもの |
| pair_swp_intrin | 4.282396 | 上記を4倍展開してSIMD化したもの |
| sorted | 7.040513 | i粒子でソートしたもの |
| sorted_intrin | 3.889211 | i粒子でソートしたものをSIMD化したもの |
| sorted_swp | 5.806741 | i粒子でソートしたものをSWPしたもの |
| sorted_swp_intrin | 3.272044 | i粒子でソートしたものをSWPしたものをSIMD化したもの |
| sorted_swp_intrin_mat_transp | 2.980598 | 上記の転置順序を変えたもの |
10%弱高速化されましたな。
コードは以下に。
まとめ
もともと、4倍展開する前にSIMD化してたものを4倍展開してから馬鹿SIMD化したので、いろいろ計算した後、最後に転置する形しか考えつかなかった。でも、プルリクもらってよく考えたら、先に転置しちゃった方が計算量少ないから速い。うーん、まだまだ修行が足りませんな。
この様子じゃ、データ構造とか工夫したり、ループのほどきかた工夫したりしたら、まだまだ高速化の余地がありそうですね。
オレはようやくのぼりはじめたばかりだからな このはてしなく遠いSIMD坂をよ・・・
こんどこそおしまい？