第4回 自動ベクトル化はどんな時に行われるか
今回の記事では、インテル® コンパイラーがどのような C/C++ ソースを自動ベクトル化できるか、いくつかのコードパターンを示して紹介します。典型的なメモリーアクセスを伴うループ
インテル® コンパイラーによる自動ベクトル化は、連続したメモリーアクセスを伴う演算が含まれるループで最も効果を発揮することは、これまでに紹介しました。次のようなコードがその代表的な例です。
for(i=0; i<MAX; i++){
A[i] += B[i] * C[i];
}
しかし、上記のループでループインデックスの増分が、+1 ではなく +2 のようにメモリーアクセスをスキップするアクセスの場合、「効率が悪いのでベクトル化しません」という様なレポートが得られるでしょう。
ループ中に関数呼び出しがある場合、原則としてベクトル化できません。しかし、呼び出し先の関数をインライン展開し、展開したコードに依存性が無いなどベクトル化の要件を満たしていればそのループはベクトル化可能です。この場合、 /Qipo (-ipo) や /Qip (-ip) オプションを指定してコンパイルします。ただし、インテル® コンパイラーに含まれる三角関数などいくつかの算術関数は、SVML（ショート・ベクター・マス・ライブラリー）でベクトル化されたバージョンが用意されています。
for(i=0; i<MAX; i++){
A[i] += func(B[i]);
}
SVML を利用する場合、mathimf.h をインクルードしてください。このファイル中にプロトタイプが宣言されています。サポートされる関数の種類については、コンパイラーのユーザーリファレンス・マニュアルの「コンパイラー・リファレンス」にある「インテルの算術ライブラリー」の項を参照してください。
また、演算を伴わない単純なメモリー操作のみを行うコードにもベクトル化を適用できます。C 標準ライブラリーで提供される strncpy() に相当する次のコードを考えてみましょう。SSE 命令は 128 ビットのレジスター（AVX では 256 ビット）へのロードとストア命令を備えており、バイト単位でメモリーを操作するよりもはるかに高速にコピーできます。
void mcopy(char *cp1, char *cp2, int n){
int i;
for(i=0; i<n; i++){
*(cp1 + i) = *(cp2 + i);
}
}
コンパイラーがこのコードをベクトル化するにあたり、不明な点が 2 つあります。 1 つは、メモリーポインター cp1 と cp2 のエリアシングの問題、そしてもう 1 つはコピーするバイト数 n が判らないことです。しかし、インテル® コンパイラーはこのような小さなバイト処理に特殊な解釈を適用します。このコードを 32 ビット版インテル® コンパイラー V12で、ベクトル化レポートオプションを指定してコンパイルしてみます。
＞ icl sample.c /Qvec-report2C:\sample.c(3): (col. 2) remark: ループがベクトル化されました。
C:\sample.c(3): (col. 2) remark: ループはスキップされました: 複数のバージョンがあります。
コンパイラーはこのコードをベクトル化しましたが、2 番目のメッセージに注目してください。コンパイラーは、この関数の2つのポインターにエリアスが無いと仮定したベクトル化コードと、エリアスがあるかループ回数がSSE命令を利用するのに適切ではない（この場合16以下）スカラーコードの複数の実行パスを生成します。その分コードサイズは少し大きくなってしまいますが、関数がどのような呼ばれ方をしても適切に処理できます。
もちろん第3回で紹介した、エリアスが無いことをコンパイラーに知らせるような機能を利用すると、その分コードパスは少なくなるので、オブジェクトサイズは小さくなります。
演算主体のループ
メモリーアクセスを伴わない演算主体のループでもベクトル化は有効です。次のコードは、数値積分によって円周率を求めるループですが、この様な同じデータ型に対し同じ演算を繰り返すループも自動ベクトル化が適用されます。この例のループの場合、ループインデックス i は i+=2 に変更され、ループ内で SIMD 命令によって 2 つの倍精度浮動小数点要素が演算されます。
for (i = 0; i < STEP; i++){
double x = (i + .5) * STEP;
sum += 4.0 / (1. + x * x);
}
ベクトル化が適用されると、ループ回数と演算数がほぼ半分になり、スループットは 2 倍になることが期待されます。インテル® Core™ i7 965 (3.2GHz/L3 8MB/HT OFF)+12GB メモリー上で動作するWindows 7 Ultimate 64bit を搭載するシステムで、ベクトル化されたコードと、スカラーコードの実行時間を計測してみたところ、それぞれ3.29 秒と6.59 秒となりました。スカラーコードは、/Qvec- オプションを指定しベクトル化を無効にしています。/Qvec- (Windows)、-no-vec (linux, Mac OS X) は、SIMD 命令を利用したスカラーコードの生成を強制します。
シーケンシャルコード
ループを伴わないコードシーケンスでもベクトル化が適用されることがあります。例えばループを開いてシーケンシャルにコードを書いたようなケースです。このような場合、ベクトル化レポートでは、「ループがベクトル化された」というメッセージの変わりに、「BLOCK がベクトル化された」というメッセージが出力されます。
次の例の各配列要素が float 型である場合、1つの SIMD 命令で処理できるデータ要素が4個分揃っているのでベクトル化されます。
A[0] = B[0] * A[0]; // 13行目 A[1] = B[1] * A[1]; A[2] = B[2] * A[2]; A[3] = B[3] * A[3];
このソースのベクトル化レポートを生成すると、次のようなメッセージが出力されます。
C:\sample.c(13): (col. 2) remark: BLOCK がベクトル化されました。最初の3行までしかない場合、データ要素が不足しているためベクトル化されません。この場合、/Qvec-report2 としてベクトル化できなかったコード領域を検出しようとしてもレポートは生成されません。ベクトル化レポートはループに対して適用されると考えてください。
ループのベクトル化を妨げる問題
これまでインテル® コンパイラーが C/C++ のループをベクトル化できるパターンと、できないパターンをいくつか紹介しましたが、ここでまとめてみます。
① 依存性があってはいけません。
for (j=1; j<MAX; j++) {
A[j] = A[j-1] + 1;
}
ループの反復間で利用するデータに依存性があってはいけません。上記の例の場合、次の関係が成り立ち、現在のループの書き込みが終了しないと、次のループは読み込みを開始できません。
A[1] = A[0] + 1 1回目のループ A[2] = A[1] + 1 2回目のループ A[3] = A[2] + 1 3回目のループ A[4] = A[3] + 1 4回目のループ … A[n] = A[n-1] + 1 n回目のループ
SIMD 命令によるベクトル化では、上記の配列 A が32ビット整数であるとすると 4 つの要素を同時に処理しようとします。ループ反復間で書き込み、参照するデータに依存関係があるとベクトル化できません。
この様な場合、ループ反復間に依存関係が無くなるようコードを書き換えることでベクトル化が可能になります。前述のループは、A[0] + 2 == (A[0] + 1 + 1) == A[1] + 1 の関係が成り立つため、次のように依存性のないループに書き換えることができます。
for (j=1; j<MAX; j++) {
A[j] = A[0] + j;
}
また、次のように書き込み先と読み込み元のエイリアスが不明な場合も依存関係がある可能性があります。このケースでは、冒頭の文字列コピーの例で紹介したように、複数のバージョンのコードが生成されます。
void func(float *z, float *x) {
int i, A = x[0];
for (i = 0; i < 100; i++)
z[i] = A * x[i];
}
すべての依存性をこのように排除できるわけではありませんが、依存性の排除はベクトル化だけではなく、他の最適化や並列化においても重要です。
② ユニットストライドではないメモリーアクセスはベクトル化できないユニットストライドではないメモリーアクセスとは、冒頭の典型的なメモリーアクセスを伴うループの例で紹介したような、非連続なメモリーアクセスのことです。ループインデックスのスキップのほか、複数次元配列の外側の次元のアクセスや、条件式を伴うメモリーアクセスがあります。
for (i=0; i<=MAX; i++) // 6行目
for (j=0; j<=MAX; j++) { // 7行目
c[i][j] += 1; // ユニットストライド 8行目
c[j][i] += 1; // 非ユニットストライド 9行目
if(a[MAX - j] == 1) // 条件式 10行目
c[j][i] = j; // 非ユニットストライド 11行目
}
8行目と9行目、そして10－11行目はそれぞれ単独行でコンパイルしてください。
上記の例で、8行目の c[i][j] はユニットストライドなアクセスですが、9行目の c[j][i] は内側のループでスキップしたメモリーアクセスとなります。バージョン 11 以降のインテル® コンパイラーでは、9行目のようなメモリーアクセスに対し自動的にループ変換を行い、ベクトル化レポートでは次のようなメッセージを生成します。
C:\sample.c(6): (col. 5) remark: PERMUTED LOOP がベクトル化されました。
この場合、6行目と7行目のループが入れ替えられ、c[j][i] のインデックス i が連続したアクセスに置き換えられます。PERMUTED LOOP とは、ループ変換されたループのことを指します。
また、10行目と11行目のようにメモリーアクセスに条件式を伴う場合、連続したアクセスとなる保証がないのでコンパイラーはベクトル化しません。レポートでは次のようなメッセージが表示されるでしょう。
C:\sample.c(7): (col. 2) remark: ループはベクトル化されませんでした: ベクトル化は可能ですが非効率です。
③ ループカウントが不明なループコンパイル時にループカウントが不明なループや実行時にカウントを判定できない場合、非標準なループとして認識されベクトル化できません。多くの場合ループカウントが不明なループに対し、最新のコンパイラーは代替えループを生成し、充分に大きなループカウントではベクトル化し、それ以外はスカラー処理を行うコードを生成します。
struct _x { int data; int bound; };
func(int *a, struct _x *x) { // 2行目
for (int i = 0; i < x->bound; i++)
a[i] = 0;
}
しかし、上記のようにループカウントに構造体要素を直接参照しているような場合、「非標準なループ」と認識し、ベクトル化できません。レポートでは次のようなメッセージが出力されます。
C:\sample.c(2): (col. 1) remark: ループはベクトル化されませんでした: 非標準のループはベクトル化候補ではありません。
このような場合、いったん整数変数にループカウントを代入し、その変数をループ式中で利用することで、コンパイラーは代替えループを生成できるようになります。
struct _x { int data; int bound; };
func(int *a, struct _x *x) {
int count = x ->bound ; // 変数に代入
for (int i = 0; i < count; i++) // 代入した変数を利用
a[i] = 0;
}
インテル® コンパイラーは、ほぼ半年ごとにバージョンアップを行っており、1年ごとにメジャーバージョンが上がっています。11.1から12.0のようにメジャーバージョンが上がる際に、大幅な機能拡張が行われており、その際にベクトル化能力も改善されてきました。2005年6月に公開されたバージョン9.0では、ベクトル化できなかった多くのケースが最新のバージョン12.0では対応できるようになっています。コンパイラーのバージョンを上げた場合、ベクトル化サポートを比較してみるといいでしょう。
さて、次回はコンパイラーにベクトル化を任せるだけではなく、明示的にプログラマーがベクトル化を記述する方法を紹介します。
