approximate double precision division using SSE2

by syoyo

SSE 命令には rcppd なる、double の逆数を近似で求める命令はありません(少なくとも SSE4 までには).
float の場合には rcpps 命令が存在し、これは 1.0 / a の近似を高速に求めることができます.

なぜ rcppd が無いのか?

rcp では、内部では除算テーブルを引くだけの実装になっているようです.
そして、実際のところ、逆数を求めるにはそれほど大きなテーブルが必要無いことが知られています.
(詳細は「ディジタル数値演算回路の実用設計」あたりを参照してください)

なのできっと、double 用にテーブルを用意するよりも、
「一旦 double を float に変換してから rcp を呼び、その結果をまた float から double に変換すればいいんじゃね?」
という考えなのかもしれません.

実際 rcppd2ps, rcpps などで検索したら、上記を行うなんともビンゴなコードが!

http://yoffy.dyndns.org/trac/yofcpp/changeset/518/trunk/yoffy/simd_type/sse2_double2.hpp

Newton 法込みで、44bits ほどの精度だそう.
というわけで、これをベースにして、divpd と cvt/rcpps + Newton 法の比較をとってみました.

ベンチマークソースコードはこちら.
http://gist.github.com/163957

Core2(Merom) 2.16GHz の結果では以下のようになりました.


Macintosh-3: $ gcc -O3 -msse2 -mfpmath=sse main.c
Macintosh-3: $ ./a.out
[SSE approx div] 0.161568 secs
[Scalar div] 0.820053 secs
SSE ver is 5.075590 times faster.

cvt/rcp で近似除算を求めるほうが、divpd で除算を求めるよりも 5 倍早いという結果になりました.
Nehalem だと、divpd, divsd のレイテンシ/スループットが向上しているので、
これよりも差は縮まるかもしれません.

Advertisements