Visual* C++、GCC*、インテル® C++ コンパイラーなど、広く利用されているコンパイラーには様々な最適化オプションが用意されています。SIMD 命令にはどれほどの効果があるのか、そして通常の IA 命令との違いは？ この回では、SIMD 命令の効果について説明します。
第2回 SIMD 命令と伝統的な IA 命令
SIMD 命令を利用すると劇的にアプリケーションの性能を改善できることがあります。Visual* C++、GCC*、インテル® C++ コンパイラーなど主要なコンパイラーは、SIMD 命令を利用した実行コードを生成することができますが、コンパイラーによってその最適化能力が異なります。利用するコンパイラーの特性を知っておくと良いでしょう。
for(int i; i<MAX; i++)
A[i] += B[i] * C[i]; // A, B, C の配列が単精度浮動小数点の場合
上記のように連続する配列データを繰り返しアクセスして演算するようなループが、最も SIMD 命令に適しています。この例で扱う配列データの型は単精度浮動小数点であり、コンパイラーは次のような実行コードを生成することが考えられます。
- 伝統的な x87 浮動小数点命令を利用した演算
- SIMD 命令を利用したスカラー演算
- SIMD 命令を利用したベクトル演算
伝統的な x87 浮動小数点命令を利用：
IA32 環境では長い間浮動小数点演算に x87 と呼ばれる命令が利用されていました。x87 命令はダイレクトにアクセスできる専用レジスターを持たず、8つのスタックを利用して演算を行います。そのため、スタック操作による制約を受ける傾向があり、一般的なループのアンロールが非常に困難です。当然1つの x87 命令は 1 つのデータしか処理できないため、上記のソースの場合、MAX 回ループを繰り返すことになります。
次のコードは、32 ビット版 Visual C++ バージョン 16 が生成したアセンブリーコードです。x87 アセンブリー命令は、先頭が “f” で始まり、コード中の fmul で 配列BとCを乗算し、fadd でその加算を行い、fstp で結果を配列 A に書き込んでいます。
$LN5@main:
mov eax, DWORD PTR _i$[ebp] // i をメモリーから読み込み
add eax, 1 // i をインクリメント
mov DWORD PTR _i$[ebp], eax // i をメモリーに退避
$LN6@main:
cmp DWORD PTR _i$[ebp], 1024 // i <=1024 を判定
jge SHORT $LN4@main // 上記が真ならループを終了
mov ecx, DWORD PTR _i$[ebp] // インデックス iを読み込み
fld DWORD PTR _B[ecx*4] // B[i] を読み込み
mov edx, DWORD PTR _i$[ebp] // インデックス iを読み込み
fmul DWORD PTR _C[edx*4] // B[i] * C[i]
mov eax, DWORD PTR _i$[ebp] // インデックス iを読み込み
fadd DWORD PTR _A[eax*4] // A += (B * C)
mov ecx, DWORD PTR _i$[ebp] // インデックス iを読み込み
fstp DWORD PTR _A[ecx*4] // A をメモリーにセーブ
jmp SHORT $LN5@main // 繰り返し
$LN4@main:
リスト1： cl sample.c /Fa で生成したアセンブリーコード
SIMD 命令を利用してスカラー演算を行う：
SSE2 などの SIMD 命令には 2 つの演算タイプがあります。XMM と呼ばれる 128 ビットレジスターには、複数のデータを格納できます。単精度浮動小数点データであれば、4 つのデータを格納することができます。データ要素を複数パックして一回の演算で処理することをベクトル演算、128 ビットレジスターの1つのデータ要素しか利用しない場合をスカラー演算と呼びます。当然ベクトル演算を行ったほうが効率良いわけです。
注：Windows 環境では、float と double のデータ精度は仮数53ビット（指数部を含めて 64 ビット）として演算してしまいます（メモリー上には、float は仮数部 23 ビットとして保存される）。float を 23 ビット単精度としてコンパイルしたい場合は、/fp:fast オプションを利用します。
次のリストは、32 ビット版 Visual C(cl)の/arch:SSE2 と /fp:fast オプションを指定して SSE2 命令を利用したアセンブリーコードです。
実際に A[i] += B[i] * C[i] を演算している命令は、mulss （乗算）と addss （加算）です。 mulss は、 MULtiply Scalar Single （単精度浮動小数点のスカラー乗算）で、addss は ADD Scalar Single（単精度浮動小数点のスカラー加算）です。つまり、128 ビットレジスターの下位 32 ビットしか使わないスカラー演算が行われています。1 命令で 1 つのデータしか処理しないわけですから x87 と効率は同じですが、SIMD 命令はダイレクトアクセスできる XMM レジスターが 8 個（64 ビットモードでは 16 個）あり、並列に実行できる内部ユニットがあるので、x87 よりは高速です。
$LN5@main:
mov eax, DWORD PTR _i$[ebp] // i をメモリーから読み込み
add eax, 1 // i をインクリメント
mov DWORD PTR _i$[ebp], eax // i をメモリーに退避
$LN6@main:
cmp DWORD PTR _i$[ebp], 1024 // i <=1024 を判定
jge SHORT $LN4@main // 上記が真ならループを終了
mov ecx, DWORD PTR _i$[ebp] // インデックス i を ecx へ設定
mov edx, DWORD PTR _i$[ebp] // インデックス i を edx へ設定
movss xmm0, DWORD PTR _B[ecx*4] // Bの1つの要素を読み込み
mulss xmm0, DWORD PTR _C[edx*4] // Cの1つの要素を読み込み
// xmm0 = B * C
mov eax, DWORD PTR _i$[ebp] // インデックス i を eax へ設定
adds xmm0, DWORD PTR _A[eax*4] // Aの1つの要素を読み込み
// xmm0 = xmm0 + A
mov ecx, DWORD PTR _i$[ebp] // インデックス i を ecx へ設定
movss DWORD PTR _A[ecx*4], xmm0 // A の要素（単精度）をメモリーへセーブ
jmp SHORT $LN5@main // ループの先頭へ
$LN4@main: // ループの出口
リスト2： cl sample.c /Fa /arch:SSE2 /fp:fast で生成したアセンブリーコード
64 ビット環境向けのコンパイラーは、浮動小数点演算をデフォルトで SSE2 命令を利用します（64 ビットをサポートする IA プロセッサーは、必ず SSE2 命令を備えている）。そのため 32 ビットのビルド環境を移行して 64 ビットコンパイラーでコンパイルすると、何もしなくても SIMD 命令が利用されることになります。
SIMD 命令を利用してベクトル演算を行う：
SIMD 命令を利用したベクトル化は命令レベルの並列性を高めますが、コンパイラーの最適化能力が影響します。コンパイラーがソースを解析しベクトル化できないと判断すると、前述の SIMD 命令を利用したスカラー演算を行うバイナリーを、もしくは SIMD 命令では効率が悪いと判断すると x87 命令を出力します。
次のアセンブリーは、32 ビット版インテル® C++ コンパイラー バージョン 12(icl)のデフォルトオプションで生成したコードです。このコードでは、ループがアンロールされている事がわかります。movaps は 16 バイト境界のメモリー領域から、16 バイトのデータ（この場合 4 つの単精度浮動小数点データ）を、128 ビットレジスターに読み込みます。
mulps は、MULtiply Packed Single（単精度浮動小数点のベクトル乗算）、addps は、ADD Packed Single（単精度浮動小数点のベクトル加算）です。これらの命令では 1 命令で、4 つのデータ要素を演算でき、ループがアンロールされているため、1 回のループで 8 個の単精度データを演算しています。
プロセッサー内部の演算動作についてはここでは触れませんが、ベクトル化されたコードでは前述のコードに比べかなり効率よく処理を行っています。コードがベクトル化されているかどうかを判断するには、アセンブリーコードを出力して SIMD 命令のサフィックスに mulps や addps のように、”px” （p はパックドの意味。x はデータ型により異なる）がついているかを指標とすると良いでしょう。
.B1.2: // ループの先頭
movaps xmm0, XMMWORD PTR [_B+edx*4] // Bを4要素読み込み
movaps xmm1, XMMWORD PTR [_B+16+edx*4] // Bの次の4要素読み込み
mulps xmm0, XMMWORD PTR [_C+edx*4] // xmm0= B * C
mulps xmm1, XMMWORD PTR [_C+16+edx*4] // xmm1= B1 * C1
addps xmm0, XMMWORD PTR [_A+edx*4] // xmm0=xmm0+A
addps xmm1, XMMWORD PTR [_A+16+edx*4] // xmm1=xmm1+A1
movaps XMMWORD PTR [_A+edx*4], xmm0 // A の4要素をメモリーへ
movaps XMMWORD PTR [_A+16+edx*4], xmm1 // 次の4要素をメモリーへ
add edx, 8 // インデックス i を 8 インクリメント
cmp edx, 1024 // i < 1024 を判定
jb .B1.2 // 真ならばループの先頭へ戻る
リスト3： icl sample.c /Faで生成したアセンブリーコード
ちなみに次のリストは、同じソースをインテルコンパイラーの AVX オプション（/QxAVX）で生成した例です。AVX では YMM という 256 ビットレジスターが利用できるため、単精度浮動小数点であれば、1 命令で 8 つのデータ要素を扱うことができます。
.B1.2:
vmovups ymm0, YMMWORD PTR [_B+edx*4] // Bを8要素読み込み
vmovups ymm3, YMMWORD PTR [_B+32+edx*4] // Bの次の8要素読み込み
vmulps ymm1, ymm0, YMMWORD PTR [_C+edx*4] // ymm1 = B * C
vmulps ymm4, ymm3, YMMWORD PTR [_C+32+edx*4] // ymm4 = B1 * C1
vaddps ymm2, ymm1, YMMWORD PTR [_A+edx*4] // ymm2 = ymm1 + A
vaddps ymm5, ymm4, YMMWORD PTR [_A+32+edx*4] // ymm5 = ymm4 + A1
vmovups YMMWORD PTR [_A+edx*4], ymm2 // A の8要素をメモリー
vmovups YMMWORD PTR [_A+32+edx*4], ymm5 // 次の8要素をメモリー
add edx, 16 // インデックス i を 16 インクリメント
cmp edx, 1024 // i < 1024 を判定
jb .B1.2 // 真ならばループの先頭へ戻る
リスト4： icl sample.c /Fa /fp:fast /QxAVX で生成したアセンブリーコード
各コンパイラーの 32 ビット環境と 64 ビット環境の SIMD デフォルトは以下のようになります；
| インテル C++ 32bit/64bit | Visual C 32bit/64bit | GCC 32bit/64bit | |
|---|---|---|---|
| 利用する命令 | IA32+SSE2/IA32+SSE2 | IA32/IA32+SSE2 | IA32/IA32+SSE2 |
また、各コンパイラーの SIMD 命令生成オプションは次のようになります；
インテル C++ 11.x と 12.0:/Qx[SSE2, SSE3, SSSE3, SSE3_ATOM,SSE4.1, SSE4.2, AVX] (Windows)
-x[SSE2, SSE3, SSSE3, SSE3_ATOM,SSE4.1, SSE4.2, AVX] (Linux & Mac OS X)
Visual C++ 15 と 16:
/arch:[SSE, SSE2, AVX] AVX は64ビット版のVC16のみ
GCC 3.x と 4.x:
-m[SSE2, SSE3, 3DNOW, SSE4a, SSE4a, SSE4.1, SSE4.2, AVX] 利用できる命令は gcc のバージョンで異なります
-mfpmath=SSE 32 ビットコンパイラーでは、このオプションも必要
特に 32 ビット環境向けのアプリケーションを開発する際には、利用するコンパイラーの SIMD オプションを有効にして評価してみてください。この例では浮動小数点データを扱いましたが、SIMD 命令は整数演算や論理演算でも利用できるので、浮動小数点演算以外の性能を改善できる場合もあります。
各コンパイラーともデフォルト以外の最適化オプションを利用すると、それぞれ最適化を凝らしたコードを生成するので、この例で紹介した例とまったく同じリストが生成されるとか限りません。次回からは、実際にいろいろなコードを利用して SIMD 命令を利用したベクトル化の可能性を考えてみましょう。
