SIMD演算
C言語でSSEやSSE2を使う方法について。ポイントがいくつかあります。
- xmmintrin.hやemmintrin.hをインクルードする
SSE命令のみならxmmintrin.h、SSE2命令も使うならemmintrin.h
- SIMD命令でアクセスするメモリは16バイトアラインメントにする
Pentium4以降に限定されますが、SSEやSSE2を使うのが良いでしょう。
これらを踏まえて、画像の各ピクセルのRGB値を反転させる例を示します。
各ピクセルは8bitのRGB値が並んでいるものとします。8bitの整数演算を行うので、8bit計算を16個並列で行うSSE2命令を使用します。
最初はSSE2のインクルードファイル。
#include <emmintrin.h>
つづいて16バイトアラインメントれたメモリを確保するために、_mm_mallocを使ってメモリを確保します。
data = _mm_malloc(SIZE,sizeof(__m128i));
__m128iは128bit(=16バイト)のデータをあらわす型です。
ピクセルの反転は255-R, 255-G,255-Bで行います。SSE2なので8bitの16個並列引き算命令が使えます。n1を255、*pInは各ピクセルのR,G,B値とすると以下のようになります。
_mm_subs_epu8 (n1 , *pIn)
コンパイルはgccの場合 -msse2オプションをつけます。
gcc -O2 -msse2 img_inv.c
以下のソースでそれぞれの処理を100回ループした結果は以下のとおり。SSEが速い。でも、16個同時だからもっと速くなると思ったけどこのくらいなのね。
通常 2.0[秒] SSE2使用 1.5[秒]
以下はソースコード。http://www.codeproject.com/KB/recipes/mmxintro.aspxを参考にしました。
/*
gcc -O2 -msse2 img_inv.c
*/
#include <stdio.h>
#include <time.h>
#include <emmintrin.h>
void InvertImage(unsigned char* pSource,
unsigned char* pDest,
int nNumberOfPixels)
{
int i;
for(i = 0; i < nNumberOfPixels; i++ ) {
*pDest = 255 - *pSource;
pDest++;
pSource++;
}
}
void InvertImage_SSE2(unsigned char* pSource,
unsigned char* pDest,
int nNumberOfPixels)
{
int nLoop = nNumberOfPixels/sizeof(__m128i);
int i;
/* ポインタとして扱う */
__m128i* pIn = (__m128i*) pSource;
__m128i* pOut = (__m128i*) pDest;
/* r0-r16のレジスタに255を設定 */
__m128i n1 = _mm_set1_epi8 (255);
for(i = 0; i < nLoop; i++ ) {
*pOut = _mm_subs_epu8 (n1 , *pIn); /* 255 - *pIn*/
/*次の16バイト*/
pIn++;
pOut++;
}
}
long diffMillis(struct timespec* t1, struct timespec* t2){
long diff = 0;
diff = (t1->tv_sec - t2->tv_sec) * 1000;
diff = diff + (t1->tv_nsec - t2->tv_nsec) / 1000000;
return diff;
}
#define SIZE 3*3200*2400
int main(void){
struct timespec t1, t2;
unsigned char *data;
unsigned char *out;
/* 16バイトアラインメントでメモリを確保 */
data = _mm_malloc(SIZE,sizeof(__m128i));
out = _mm_malloc(SIZE,sizeof(__m128i));
clock_gettime(CLOCK_REALTIME, &t1);
InvertImage(data,out,SIZE);
clock_gettime(CLOCK_REALTIME, &t2);
printf("InvertImage=%d\n",diffMillis(&t2,&t1));
clock_gettime(CLOCK_REALTIME, &t1);
InvertImage_SSE2(data,out,SIZE);
clock_gettime(CLOCK_REALTIME, &t2);
printf("InvertImage_SSE2=%d\n",diffMillis(&t2,&t1));
}アセンブラコード
gcc -Sでアセンブラのコードを見てみました。
まずは、通常計算。あれ？notbで計算してる！反転命令使ってるじゃん。恐るべし最適化。
L5: movzbl (%ebx), %eax addl $1, %ebx notb %al movb %al, (%ecx) addl $1, %ecx subl $1, %edx jne L5
つづいてSSE2。こちらはコードに書いたとおり_mm_subs_epu8がpsubusbとなっている。通常の演算がnotbを使ってたんで、XORを計算する_mm_xor_si128を使ってみたけど、実行速度は変わらず。
L15: movdqa %xmm1, %xmm0 psubusb (%ebx), %xmm0 addl $16, %ebx movdqa %xmm0, (%ecx) addl $16, %ecx subl $1, %eax jne L15
SSE2使用で1.5秒はあんまり速くない気がする。このデータ量だと512Mピクセル/sの計算を行っていることになる。メモリはread/write行うので 512×3×2=3G/s の帯域を使っている。このくらいが限界なのか？
CPUキャッシュに入るくらいのデータ量にすると、SSE2使用が3Gピクセル/s、通常計算が390Mピクセル/sとなり7倍程度の差が出た。メモリ帯域がネックのようだ。
ちょっと思い当たって以下のところを見直してみた。
*pOut = _mm_subs_epu8 (n1 , *pIn); /* 255 - *pIn*/
これ、元データを破壊していいなら以下のように書ける
*pIn = _mm_subs_epu8 (n1 , *pIn); /* 255 - *pIn*/
こっちの方がキャッシュが効くので速くないか？と思って試したら、かなり速くなった。
通常 2.0[秒] SSE2使用 1.2[秒]
もうちょっと調べたらSSE2にはキャッシュにデータを残さないでメモリに転送する命令があった。これを使ったコードは以下のとおり。
n2 = _mm_subs_epu8(n1, *pIn); /* 255 - *pIn*/ _mm_stream_si128(pIn,n2);
結果としては 1.2[秒]→0.9[秒]に！ SSEは奥が深いな。
