昔話
それは昔々のこと。
x86には浮動小数点演算を行う手段がなく、外付けの浮動小数点演算ユニットを接続するという手法で、浮動小数点演算を実現していたのであった。
x87と呼ばれたそれはとてもエクセレントなシステムで…という話はwikipediaに譲ろう。
https://ja.wikipedia.org/wiki/Intel_8087
重要なのは、x87が内部表現として80bitの拡張倍精度を使っている、ということ。
これのおかげで、x87においては、確かに~~doubleのほうが速かった
(floatだとdoubleへのキャストコストが発生するため)~~
嘘だろそれ。ASM見たら別にキャストとかしてなかったわ。
どっちかというと丸めの影響で精度が異なることのほうが重要だわ。
改めて調べてみると、doubleが速いとされている資料についてはあんまりないことに気付く。
(同等としている資料はherumiさんの http://homepage1.nifty.com/herumi/prog/prog90.html がある)
速いってのはそれがデマだったのでは…？
(2016/06/13追記: コメントにてご指摘いただきましたが、確かにfloatがdoubleより遅い条件があります。詳しくは後述)
しかしx87の登場からすでに30年(!)
x86も時代が移り変わり、x87は内部に統合され、MMX(こいつももはや死語),SSE,AVXと色々な命令が増えてきた。
この時代にあって、doubleがfloatより本当に速いのかどうか、調べてみよう。
ソースコード
gistに置いてあるので適当にとっていってほしい。
https://gist.github.com/telmin/846e1ee9ecddd4fcc063
仮説
速くないと一刀両断するのは簡単なので、doubleとfloatでどういう性能差が発生するかまず考えてみる。
キャッシュのはなし
floatとdoubleはご存知のとおり32bitと64bitで、使用するデータ量は倍になる。
とすると、当然キャッシュヒットしないような状態で計算するとあっさり性能は倍だろう。これはあんまり面白くないのと、いろいろと他でややこしくなるので、全データがキャッシュに乗るような状態を考えよう。
SSE命令のはなし
現在、浮動小数点SIMD命令として、x86(_64)には、SSE, AVXが用意されている。
話を簡単にするためにSSEのみにフォーカスする。
最近のAuto-Vectorizationは結構優秀で、簡単な演算は割と勝手にベクタライズしてくれる。ありがたいことだ。それでもまだ手で書きたい。一日一回感謝のIntrinsics!!!
でまぁ、floatとdoubleのmulでどのようにレイテンシとスループットが違うかというと
Intel Intrinsics guide
https://software.intel.com/sites/landingpage/IntrinsicsGuide/
より
- __m128 _mm_mul_ps (__m128 a, __m128 b)
| Architecture | Latency | Throughput |
|---|---|---|
| Haswell | 5 | 0.5 |
| Ivy Bridge | 5 | 1 |
| Sandy Bridge | 5 | 1 |
| Westmere | 4 | 1 |
| Nehalem | 4 | 1 |
- __m128d _mm_mul_pd (__m128d a, __m128d b)
| Architecture | Latency | Throughput |
|---|---|---|
| Haswell | 5 | 0.5 |
| Ivy Bridge | 5 | 1 |
| Sandy Bridge | 5 | 1 |
| Westmere | 4 | 1 |
| Nehalem | 4 | 1 |
となり、全く一緒になる。
なのであとは純粋にbit長に依存して、4Flop, 2Flopとなるわけだ。
あと、SSEにはxmmレジスタを入力としながらも先頭32bit, 64bitのみを演算として使用する命令もあり、普段x86_64で浮動小数点を使うとこちらが使用される。あとベクタライズしろって言ってんのにこいつを吐くこともある。
これが、それぞれ
- __m128 _mm_mul_ss (__m128 a, __m128 b)
- __m128 _mm_mul_sd (__m128 a, __m128 b)
レイテンシ、スループットはps, pdと同じ(当たり前)
x87命令のはなし
SSEのときと同様に、x87命令のレイテンシ、スループットを考えてみる。
x87の浮動小数点乗算命令はfmulである。
レイテンシ、スループットはhttp://www.agner.org/optimize/instruction_tables.pdf より
| Architecture | Latency | Throughput |
|---|---|---|
| Haswell | 5 | 1 |
| Ivy Bridge | 5 | 1 |
| Sandy Bridge | 5 | 1 |
となっている。
なんか全体的に同じ感じするわな。
つまり
ベクタライズしていればfloatとdoubleが約2倍の性能になるはず。
してない場合はx87、SSEのいずれにおいても、doubleとfloatは同等の性能になるはず。
x87の場合は、float -> doubleへのキャストコストが発生するからdoubleのほうが若干速いはず 嘘。同じ速度ぐらいになるはず。
実験環境 & 実験条件
環境
MacbookPro Retina 15(mid 2012)
CPU : Corei7-3720QM IvyBridge
RAM : 8GB
コンパイラ : Apple LLVM version 6.0 (clang-600.0.54) (based on LLVM 3.5svn)
条件
float, double 4096要素に対してaxpyを実施
中間変数含め、使うのは4096 * 4要素
axpyを100回実行してかかった時間を合計する。
x87命令動かすのに32bitバイナリじゃないとちゃんと動かなかったので以下はm32オプションを付与しているものとする。
つかなんで動かないのか誰か教えてプリーズ。
ABI絡みだった気がする…？
実験結果
とりあえずSSE
$ clang++ -O2 double.cpp -std=c++11 -m32
float time 185000 nsec
double time 390000 nsec
x87
$ clang++ -O2 double.cpp -std=c++11 -m32 -mno-sse
float time 561000 nsec
double time 549000 nsec
SSEでベクタライズ止めてみる
$ clang++ -O2 double.cpp -std=c++11 -m32 -fno-vectorize
float time 549000 nsec
double time 553000 nsec
表にまとめると以下のようになる。
| float | double|
----|----|----
SSE| 185000 | 390000 |
x87| 561000 | 549000 |
SSE(no-vectorize)| 549000 | 553000 |
x87でfloatが若干遅いのはなんだろうか。何度か実行してみると結果が若干揺れるからそのせいかな…
まとめ
あんまり優位な差は出なかったので結論が微妙になりそうな気がしないでもないが、少なくともSSEとかがあるこんにちでは、精度的に大丈夫な場合はfloatを使うと、SIMDの恩恵にあずかれるからそっちのほうがいいんじゃないかとか思う。
あと、x87は頑張らないと出ないので、もう拡張倍精度のことは忘れたほうがいい。マジで。
おまけ
ARMでの結果は以下
http://d.hatena.ne.jp/nakamura001/20090301/1235924387
(2009年なので、最新だとどうなってるかはわからないけどあんまり変わらんでしょ)
おまけその2
ちゃんとSSE命令とか出てるぅ？
SSE版
LBB0_17: ## %vector.body116
## Parent Loop BB0_13 Depth=1
## => This Inner Loop Header: Depth=2
movups (%eax,%ebx,4), %xmm0
movups (%ecx,%ebx,4), %xmm1
mulps %xmm0, %xmm1
movups (%edx,%ebx,4), %xmm0
addps %xmm1, %xmm0
movups %xmm0, (%esi,%ebx,4)
addl $4, %ebx
cmpl $4096, %ebx ## imm = 0x1000
jne LBB0_17
x87版
LBB0_14: ## Parent Loop BB0_13 Depth=1
## => This Inner Loop Header: Depth=2
flds (%eax,%ebx,4)
fmuls (%ecx,%ebx,4)
fadds (%edx,%ebx,4)
fstps (%esi,%ebx,4)
incl %ebx
cmpl $4096, %ebx ## imm = 0x1000
jne LBB0_14
出てるっぽい。
おまけその3
floatがdoubleより遅くなるシチュエーションは存在した…！
K&Rの頃のC言語は、floatはdoubleに昇格して演算を行うルールであった。
(fujitanozomu さんにコメントいただきました。ありがとうございました)
https://docs.oracle.com/cd/E19957-01/806-4836/compat.html
これでも記載されている通り、
| K&R | ANSI/ISO C | |
|---|---|---|
| 単精度計算と倍精度計算 | 浮動小数点式のオペランドを double に拡張する。float を返すように宣言された関数の戻り値は、常に double に拡張される。 | float の演算を単精度計算で行うことができる。このような関数に float の戻り型を使用できる。 |
となっているので、ANSIの規格化を以って、ちゃんとfloatが(言語仕様的には)使えるようになった模様。
あと、現在でも下手するとこの型昇格は踏むらしい。それが以下のような条件。
extern float fa, fb, fc, fd;
extern double da, db, dc, dd;
extern float fcomp();
extern double dcomp();
hoge()
{
fd = fcomp(fa, fb, fc);
}
piyo()
{
dd = dcomp(da, db, dc);
}
$ clang -O2 -m32 -S hoge.c -o -
_hoge: # @hoge
subl $24, %esp
movss _fa, %xmm0 # xmm0 = mem[0],zero,zero,zero
cvtss2sd %xmm0, %xmm0
movss _fb, %xmm1 # xmm1 = mem[0],zero,zero,zero
cvtss2sd %xmm1, %xmm1
movss _fc, %xmm2 # xmm2 = mem[0],zero,zero,zero
cvtss2sd %xmm2, %xmm2
movsd %xmm2, 16(%esp)
movsd %xmm1, 8(%esp)
movsd %xmm0, (%esp)
calll _fcomp
fstps _fd
addl $24, %esp
retl
_piyo: # @piyo
subl $24, %esp
movsd _da, %xmm0 # xmm0 = mem[0],zero
movsd _db, %xmm1 # xmm1 = mem[0],zero
movsd _dc, %xmm2 # xmm2 = mem[0],zero
movsd %xmm2, 16(%esp)
movsd %xmm1, 8(%esp)
movsd %xmm0, (%esp)
calll _dcomp
fstpl _dd
addl $24, %esp
retl
float fcomp(a, b, c)
float a;
float b;
float c;
{
return a + b * c;
}
double dcomp(a, b, c)
double a;
double b;
double c;
{
return a + b * c;
}
$ clang -O2 -m32 -S piyo.c -o -
_fcomp: # @fcomp
pushl %eax
movsd 24(%esp), %xmm0 # xmm0 = mem[0],zero
movsd 16(%esp), %xmm1 # xmm1 = mem[0],zero
movsd 8(%esp), %xmm2 # xmm2 = mem[0],zero
cvtsd2ss %xmm2, %xmm2
cvtsd2ss %xmm1, %xmm1
cvtsd2ss %xmm0, %xmm0
mulss %xmm1, %xmm0
addss %xmm2, %xmm0
movss %xmm0, (%esp)
flds (%esp)
popl %eax
retl
_dcomp: # @dcomp
pushl %ebp
movl %esp, %ebp
andl $-8, %esp
subl $8, %esp
movsd 16(%ebp), %xmm0 # xmm0 = mem[0],zero
mulsd 24(%ebp), %xmm0
addsd 8(%ebp), %xmm0
movsd %xmm0, (%esp)
fldl (%esp)
movl %ebp, %esp
popl %ebp
retl
現代のコンパイラでも、外部関数についてプロトタイプ宣言を行わず、昔の K&R スタイルの仮引数のない外部宣言のみで float の値を引数として与えて関数呼び出しを行うと、float の値は一旦 double に昇格されてから引数として渡されるコードが生成される
とご指摘いただいたので、この辺は気を付けないといけなさそう。