1000個の32ビット整数の配列の場合、Intel Sandybridgeのループで@hirschhornsalzのメソッドを使用して、シングルスレッドの約1.4倍の小さなスピードアップを得ることができました。intの60kiBバッファーを使用すると、スピードアップは約1.37です。8MiBのintを使用しても、スピードアップは1.13のままです。（3.8GHzターボでi5-2500k、DDR3-1600を使用。）
小さい要素（int16_tまたは、、uint8_tまたは符号なしバージョン）は、ベクトルごとの要素数が2倍になるごとに、シフト/追加の追加ステージを取ります。オーバーフローは悪いので、SSEに大きな利点があるとしても、すべての要素の合計を保持できないデータ型を使用しようとしないでください。
#include <immintrin.h>
// In-place rewrite an array of values into an array of prefix sums.
// This makes the code simpler, and minimizes cache effects.
int prefix_sum_sse(int data[], int n)
{
// const int elemsz = sizeof(data[0]);
#define elemsz sizeof(data[0]) // clang-3.5 doesn't allow compile-time-const int as an imm8 arg to intrinsics
__m128i *datavec = (__m128i*)data;
const int vec_elems = sizeof(*datavec)/elemsz;
// to use this for int8/16_t, you still need to change the add_epi32, and the shuffle
const __m128i *endp = (__m128i*) (data + n - 2*vec_elems); // don't start an iteration beyond this
__m128i carry = _mm_setzero_si128();
for(; datavec <= endp ; datavec += 2) {
IACA_START
__m128i x0 = _mm_load_si128(datavec + 0);
__m128i x1 = _mm_load_si128(datavec + 1); // unroll / pipeline by 1
// __m128i x2 = _mm_load_si128(datavec + 2);
// __m128i x3;
x0 = _mm_add_epi32(x0, _mm_slli_si128(x0, elemsz)); // for floats, use shufps not bytewise-shift
x1 = _mm_add_epi32(x1, _mm_slli_si128(x1, elemsz));
x0 = _mm_add_epi32(x0, _mm_slli_si128(x0, 2*elemsz));
x1 = _mm_add_epi32(x1, _mm_slli_si128(x1, 2*elemsz));
// more shifting if vec_elems is larger
x0 = _mm_add_epi32(x0, carry); // this has to go after the byte-shifts, to avoid double-counting the carry.
_mm_store_si128(datavec +0, x0); // store first to allow destructive shuffle (non-avx pshufb if needed)
x1 = _mm_add_epi32(_mm_shuffle_epi32(x0, _MM_SHUFFLE(3,3,3,3)), x1);
_mm_store_si128(datavec +1, x1);
carry = _mm_shuffle_epi32(x1, _MM_SHUFFLE(3,3,3,3)); // broadcast the high element for next vector
}
// FIXME: scalar loop to handle the last few elements
IACA_END
return data[n-1];
#undef elemsz
}
int prefix_sum_simple(int data[], int n)
{
int sum=0;
for (int i=0; i<n ; i++) {
IACA_START
sum += data[i];
data[i] = sum;
}
IACA_END
return sum;
}
// perl -we '$n=1000; sub rnlist($$) { return map { int rand($_[1]) } ( 1..$_[0] );} @a=rnlist($n,127); $"=", "; print "$n\n@a\n";'
int data[] = { 51, 83, 126, 11, 20, 63, 113, 102,
126,67, 83, 113, 86, 123, 30, 109,
97, 71, 109, 86, 67, 60, 47, 12,
/* ... */ };
int main(int argc, char**argv)
{
const int elemsz = sizeof(data[0]);
const int n = sizeof(data)/elemsz;
const long reps = 1000000 * 1000 / n;
if (argc >= 2 && *argv[1] == 'n') {
for (int i=0; i < reps ; i++)
prefix_sum_simple(data, n);
}else {
for (int i=0; i < reps ; i++)
prefix_sum_sse(data, n);
}
return 0;
}
リストをバイナリにコンパイルして、n=1000でテストします。（はい、実際にループしていることを確認しました。ベクトルテストまたは非ベクトルテストを無意味にするコンパイル時のショートカットは使用していません。）
3オペランドの非破壊ベクトル命令を取得するためにAVXをコンパイルすると、多くのmovdqa命令が節約されますが、節約されるサイクルはごくわずかであることに注意してください。これは、shuffleとvector-int-addの両方がSnB / IvBのポート1と5でのみ実行できるため、port0にはmov命令を実行するための十分な予備サイクルがあるためです。uop-cacheスループットのボトルネックが、AVX以外のバージョンがわずかに遅い理由である可能性があります。（これらの追加のmov命令はすべて、最大3.35 insn / cycleをプッシュします）。フロントエンドはサイクルの4.54％しかアイドル状態ではないため、ほとんど追いついていない。
gcc -funroll-loops -DIACA_MARKS_OFF -g -std=c11 -Wall -march=native -O3 prefix-sum.c -mno-avx -o prefix-sum-noavx
# gcc 4.9.2
################# SSE (non-AVX) vector version ############
$ ocperf.py stat -e task-clock,cycles,instructions,uops_issued.any,uops_dispatched.thread,uops_retired.all,uops_retired.retire_slots,stalled-cycles-frontend,stalled-cycles-backend ./prefix-sum-noavx
perf stat -e task-clock,cycles,instructions,cpu/event=0xe,umask=0x1,name=uops_issued_any/,cpu/event=0xb1,umask=0x1,name=uops_dispatched_thread/,cpu/event=0xc2,umask=0x1,name=uops_retired_all/,cpu/event=0xc2,umask=0x2,name=uops_retired_retire_slots/,stalled-cycles-frontend,stalled-cycles-backend ./prefix-sum-noavx
Performance counter stats for './prefix-sum-noavx':
206.986720 task-clock (msec) # 0.999 CPUs utilized
777,473,726 cycles # 3.756 GHz
2,604,757,487 instructions # 3.35 insns per cycle
# 0.01 stalled cycles per insn
2,579,310,493 uops_issued_any # 12461.237 M/sec
2,828,479,147 uops_dispatched_thread # 13665.027 M/sec
2,829,198,313 uops_retired_all # 13668.502 M/sec (unfused domain)
2,579,016,838 uops_retired_retire_slots # 12459.818 M/sec (fused domain)
35,298,807 stalled-cycles-frontend # 4.54% frontend cycles idle
1,224,399 stalled-cycles-backend # 0.16% backend cycles idle
0.207234316 seconds time elapsed
------------------------------------------------------------
######### AVX (same source, but built with -mavx). not AVX2 #########
$ ocperf.py stat -e task-clock,cycles,instructions,uops_issued.any,uops_dispatched.thread,uops_retired.all,uops_retired.retire_slots,stalled-cycles-frontend,stalled-cycles-backend ./prefix-sum-avx
Performance counter stats for './prefix-sum-avx':
203.429021 task-clock (msec) # 0.999 CPUs utilized
764,859,441 cycles # 3.760 GHz
2,079,716,097 instructions # 2.72 insns per cycle
# 0.12 stalled cycles per insn
2,054,334,040 uops_issued_any # 10098.530 M/sec
2,303,378,797 uops_dispatched_thread # 11322.764 M/sec
2,304,140,578 uops_retired_all # 11326.509 M/sec
2,053,968,862 uops_retired_retire_slots # 10096.735 M/sec
240,883,566 stalled-cycles-frontend # 31.49% frontend cycles idle
1,224,637 stalled-cycles-backend # 0.16% backend cycles idle
0.203732797 seconds time elapsed
------------------------------------------------------------
################## scalar version (cmdline arg) #############
$ ocperf.py stat -e task-clock,cycles,instructions,uops_issued.any,uops_dispatched.thread,uops_retired.all,uops_retired.retire_slots,stalled-cycles-frontend,stalled-cycles-backend ./prefix-sum-avx n
Performance counter stats for './prefix-sum-avx n':
287.567070 task-clock (msec) # 0.999 CPUs utilized
1,082,611,453 cycles # 3.765 GHz
2,381,840,355 instructions # 2.20 insns per cycle
# 0.20 stalled cycles per insn
2,272,652,370 uops_issued_any # 7903.034 M/sec
4,262,838,836 uops_dispatched_thread # 14823.807 M/sec
4,256,351,856 uops_retired_all # 14801.249 M/sec
2,256,150,510 uops_retired_retire_slots # 7845.650 M/sec
465,018,146 stalled-cycles-frontend # 42.95% frontend cycles idle
6,321,098 stalled-cycles-backend # 0.58% backend cycles idle
0.287901811 seconds time elapsed
------------------------------------------------------------
Haswellはほぼ同じであるはずですが、シャッフルはポート1ではなくポート5でのみ実行できるため、クロックごとに少し遅くなる可能性があります（vector-intaddはHaswellではまだp1/ 5です）。
OTOH、IACAは、コンパイルせずにコンパイルした場合、Haswellは1回の反復でSnBよりもわずかに高速になると考え-funroll-loopsています（これはSnBで役立ちます）。Haswellはport6でブランチを実行できますが、SnBブランチではport5にあり、すでに飽和状態になっています。
# compile without -DIACA_MARKS_OFF
$ iaca -64 -mark 1 -arch HSW prefix-sum-avx
Intel(R) Architecture Code Analyzer Version - 2.1
Analyzed File - prefix-sum-avx
Binary Format - 64Bit
Architecture - HSW
Analysis Type - Throughput
*******************************************************************
Intel(R) Architecture Code Analyzer Mark Number 1
*******************************************************************
Throughput Analysis Report
--------------------------
Block Throughput: 6.20 Cycles Throughput Bottleneck: Port5
Port Binding In Cycles Per Iteration:
---------------------------------------------------------------------------------------
| Port | 0 - DV | 1 | 2 - D | 3 - D | 4 | 5 | 6 | 7 |
---------------------------------------------------------------------------------------
| Cycles | 1.0 0.0 | 5.8 | 1.4 1.0 | 1.4 1.0 | 2.0 | 6.2 | 1.0 | 1.3 |
---------------------------------------------------------------------------------------
N - port number or number of cycles resource conflict caused delay, DV - Divider pipe (on port 0)
D - Data fetch pipe (on ports 2 and 3), CP - on a critical path
F - Macro Fusion with the previous instruction occurred
* - instruction micro-ops not bound to a port
^ - Micro Fusion happened
# - ESP Tracking sync uop was issued
@ - SSE instruction followed an AVX256 instruction, dozens of cycles penalty is expected
! - instruction not supported, was not accounted in Analysis
| Num Of | Ports pressure in cycles | |
| Uops | 0 - DV | 1 | 2 - D | 3 - D | 4 | 5 | 6 | 7 | |
---------------------------------------------------------------------------------
| 1 | | | 1.0 1.0 | | | | | | | vmovdqa xmm2, xmmword ptr [rax]
| 1 | 1.0 | | | | | | | | | add rax, 0x20
| 1 | | | | 1.0 1.0 | | | | | | vmovdqa xmm3, xmmword ptr [rax-0x10]
| 1 | | | | | | 1.0 | | | CP | vpslldq xmm1, xmm2, 0x4
| 1 | | 1.0 | | | | | | | | vpaddd xmm2, xmm2, xmm1
| 1 | | | | | | 1.0 | | | CP | vpslldq xmm1, xmm3, 0x4
| 1 | | 1.0 | | | | | | | | vpaddd xmm3, xmm3, xmm1
| 1 | | | | | | 1.0 | | | CP | vpslldq xmm1, xmm2, 0x8
| 1 | | 1.0 | | | | | | | | vpaddd xmm2, xmm2, xmm1
| 1 | | | | | | 1.0 | | | CP | vpslldq xmm1, xmm3, 0x8
| 1 | | 1.0 | | | | | | | | vpaddd xmm3, xmm3, xmm1
| 1 | | 0.9 | | | | 0.2 | | | CP | vpaddd xmm1, xmm2, xmm0
| 2^ | | | | | 1.0 | | | 1.0 | | vmovaps xmmword ptr [rax-0x20], xmm1
| 1 | | | | | | 1.0 | | | CP | vpshufd xmm1, xmm1, 0xff
| 1 | | 0.9 | | | | 0.1 | | | CP | vpaddd xmm0, xmm1, xmm3
| 2^ | | | 0.3 | 0.3 | 1.0 | | | 0.3 | | vmovaps xmmword ptr [rax-0x10], xmm0
| 1 | | | | | | 1.0 | | | CP | vpshufd xmm0, xmm0, 0xff
| 1 | | | | | | | 1.0 | | | cmp rax, 0x602020
| 0F | | | | | | | | | | jnz 0xffffffffffffffa3
Total Num Of Uops: 20
ところで、gccは、ループカウンターがあり、を実行している場合でも、1レジスタアドレッシングモードを使用するようにループをコンパイルしましたload(datavec + i + 1)。それが最高のコード、特にです。2レジスタアドレッシングモードがマイクロヒューズできないSnBファミリでは、clangの利益のために、ソースをそのループ条件に変更します。
