19

i7で最も効率的な方法で、浮動ベクトルとビットベクトルの間の内積を計算しようとしています。実際には、128 次元または 256 次元のベクトルでこの操作を行っていますが、説明のために、問題を説明するために 64 次元のコードを書きましょう。

// a has 64 elements. b is a bitvector of 64 dimensions.
float dot(float *restrict a, uint64_t b) {
    float sum = 0;
    for(int i=0; b && i<64; i++, b>>=1) {
        if (b & 1) sum += a[i];
    }
    return sum;
}

これはもちろん機能しますが、問題は、これがプログラム全体のタイム クリティカルな場所であるため (50 分間の実行で 95% の CPU 時間を消費する)、どうしても高速化する必要があることです。

私の推測では、上記の分岐はゲーム キラーです (順不同の実行を防ぎ、不適切な分岐予測を引き起こします)。ここでベクトル命令を使用して役立つかどうかはわかりません。-std=c99 -march=native -mtune=native -Ofast -funroll-loops で gcc 4.8 を使用すると、現在この出力が得られます

    movl    $4660, %edx
    movl    $5, %ecx
    xorps   %xmm0, %xmm0
    .p2align 4,,10
    .p2align 3
.L4:
    testb   $1, %cl
    je  .L2
    addss   (%rdx), %xmm0
.L2:
    leaq    4(%rdx), %rax
    shrq    %rcx
    testb   $1, %cl
    je  .L8
    addss   4(%rdx), %xmm0
.L8:
    shrq    %rcx
    testb   $1, %cl
    je  .L9
    addss   4(%rax), %xmm0
.L9:
    shrq    %rcx
    testb   $1, %cl
    je  .L10
    addss   8(%rax), %xmm0
.L10:
    shrq    %rcx
    testb   $1, %cl
    je  .L11
    addss   12(%rax), %xmm0
.L11:
    shrq    %rcx
    testb   $1, %cl
    je  .L12
    addss   16(%rax), %xmm0
.L12:
    shrq    %rcx
    testb   $1, %cl
    je  .L13
    addss   20(%rax), %xmm0
.L13:
    shrq    %rcx
    testb   $1, %cl
    je  .L14
    addss   24(%rax), %xmm0
.L14:
    leaq    28(%rax), %rdx
    shrq    %rcx
    cmpq    $4916, %rdx
    jne .L4
    ret

編集データを並べ替えても問題ありません (並べ替えがすべてのパラメーターで同じである限り)。順序は関係ありません。

Chris Dodd の SSE2 コードの 3 倍以上の速度で動作するものがあるかどうか疑問に思っています。

新しいメモ: AVX/AVX2 コードも歓迎です!

編集 2 ビットベクトルが与えられた場合、128 (または 256 ビットの場合は 256) の異なる float ベクトルで乗算する必要があります (したがって、一度に複数の float ベクトルを使用しても問題ありません)。これがプロセス全体です。プロセス全体をスピードアップするものも大歓迎です!

4

8 に答える 8

16

最善の策は、一度に 4 つの float を操作する SSE ps 命令を使用することです。float 0.0 がすべて 0 ビットであるという事実を利用して、andps 命令を使用して不要な要素をマスクすることができます。

#include <stdint.h>
#include <xmmintrin.h>

union {
    uint32_t i[4];
    __m128   xmm;
} mask[16] = {
 {  0,  0,  0,  0 },
 { ~0,  0,  0,  0 },
 {  0, ~0,  0,  0 },
 { ~0, ~0,  0,  0 },
 {  0,  0, ~0,  0 },
 { ~0,  0, ~0,  0 },
 {  0, ~0, ~0,  0 },
 { ~0, ~0, ~0,  0 },
 {  0,  0,  0, ~0 },
 { ~0,  0,  0, ~0 },
 {  0, ~0,  0, ~0 },
 { ~0, ~0,  0, ~0 },
 {  0,  0, ~0, ~0 },
 { ~0,  0, ~0, ~0 },
 {  0, ~0, ~0, ~0 },
 { ~0, ~0, ~0, ~0 },
};

float dot(__m128 *a, uint64_t b) {
    __m128 sum = { 0.0 };
    for (int i = 0; i < 16; i++, b>>=4)
        sum += _mm_and_ps(a[i], mask[b&0xf].xmm);
    return sum[0] + sum[1] + sum[2] + sum[3];
}

マスクに多数の 0 が含まれると予想される場合は、0 を短縮した方が速い場合があります。

for (int i = 0; b; i++, b >>= 4)
    if (b & 0xf)
        sum += _mm_and_ps(a[i], mask[b&0xf].xmm);

ただし、b がランダムな場合、これは遅くなります。

于 2013-04-17T06:50:16.620 に答える
4

dataでわずかに異なる順列を許可するかfloat data[128]、 のビットマスクで対応する順列を行う場合、__m128 mask;上記の Chris Dodd によって提案されたアルゴリズムをわずかに改善できます。(マスクの並べ替えに必要な時間をカウントしないと、この実装 (+ オーバーヘッド) は約 25% 高速になります)。もちろん、これはコメントで提供された私のアイデアの簡単なドラフトです。

union {
    unsigned int i[4];
    float f[4];
    __m128   xmm;
} mask = { 0xFF00FF00, 0xF0F0F0F0, 0xCCCCCCCC, 0xAAAAAAAA };

float dot2(__m128 *a, __m128 mask);
// 20M times 1.161s

float dotref(__m128 *a, unsigned int *mask) // 20M times 8.174s
{
    float z=0.0f;
    int i;
    for (i=0;i<32;i++) {
       if (mask[0] & (0x80000000U >> i)) z+= a[i][0];
       if (mask[1] & (0x80000000U >> i)) z+= a[i][1];
       if (mask[2] & (0x80000000U >> i)) z+= a[i][2];
       if (mask[3] & (0x80000000U >> i)) z+= a[i][3];       
    }
   return z;
}

対応するアセンブラーの実装は次のようになります。

dot2:
    // warm up stage: fill in initial data and
    // set up registers
    pxor %xmm1, %xmm1      ;; // clear partial sum1
    pxor %xmm2, %xmm2      ;; // clear partial sum2
    movaps (%rdi), %xmm3   ;; // register warm up stage1
    movaps 16(%rdi), %xmm4 ;; // next 4 values
    pxor %xmm5, %xmm5
    pxor %xmm6, %xmm6
    lea 32(%rdi), %rdi
    movl $16, %ecx            ;; // process 2x4 items per iteration (total=128)
a:  ;; // inner loop -- 2 independent data paths
    blendvps %xmm3, %xmm5
    pslld $1, %xmm0
    movaps (%rdi), %xmm3   
    blendvps %xmm4, %xmm6
    pslld $1, %xmm0
    movaps 16(%rdi), %xmm4
    addps %xmm5, %xmm1
    pxor  %xmm5, %xmm5
    addps %xmm6, %xmm2
    pxor  %xmm6, %xmm6
    lea 32(%rdi), %rdi
    loop a
 ;; // cool down stage: gather results (xmm0 = xmm1+xmm2)
 ;; // in beautiful world this stage is interleaved
 ;; // with the warm up stage of the next block
    addps %xmm2, %xmm1
    movaps  %xmm1, %xmm2
    movaps  %xmm1, %xmm0
    shufps  $85, %xmm1, %xmm2
    addss   %xmm2, %xmm0
    movaps  %xmm1, %xmm2
    unpckhps %xmm1, %xmm2
    shufps  $255, %xmm1, %xmm1
    addss   %xmm2, %xmm0
    addss   %xmm1, %xmm0
    ret
于 2013-04-19T13:13:45.033 に答える
3

試してみるべきことがいくつかあります。

CMOVブランチの代わりにコンパイラを使用するようにしてください。(このように共用体を使用することは、C11 では明確に定義されていますが、C++11 では定義されていないことに注意してください。 )

union {
    int i;
    float f;
} u;
u.i = 0;
if (b & 1) {
    u.f = a[i];
}
sum += u.f;

分岐の代わりに乗算を使用します。

sum += (b & 1) * a[i];

いくつかの合計を保持し、それらを最後に追加して、データ フローの依存関係を減らします。(上記の提案のいずれかをこれと組み合わせることができます。)

float sum0 = 0, sum1 = 0, sum2 = 0, sum3 = 0;
for (int i = 0; i < 64; i += 4; b >>= 4) {
    if (b & 1) sum0 += a[i];
    if (b & 2) sum1 += a[i+1];
    if (b & 4) sum2 += a[i+2];
    if (b & 8) sum3 += a[i+3];
}
return sum0 + sum1 + sum2 + sum3;

一度にいくつかのビットを処理して、分岐の数を減らします。

for (int i = 0; i < 64; i += 4, b >>= 4) {
    switch (b & 0xf) {
        case 0:
            break;
        case 1:
            sum += a[i];
            break;
        case 2:
            sum += a[i + 1];
            break;
        case 3:
            sum += a[i] + a[i+1];
            break;
        case 4:
            sum += a[i+2];
            break;
        // etc. for cases up to and including 15
    }
}

複数の合計を保持し、合計ごとに一度に複数のビットを処理できます。その場合、おそらくマクロまたはインライン関数を使用して、それを 4 回呼び出したいと思うでしょう。

于 2013-04-17T04:39:52.403 に答える
1

i7 を使用している場合、SSE4.1 が_mm_dp_psあり、内積を行う組み込み関数を使用できます。サンプルコードでは、次のようになります

#include <stdint.h>
#include <immintrin.h>

const float fltmask[][4] =
  {{0, 0, 0, 0}, {0, 0, 0, 1}, {0, 0, 1, 0}, {0, 0, 1, 1},
   {0, 1, 0, 0}, {0, 1, 0, 1}, {0, 1, 1, 0}, {0, 1, 1, 1},
   {1, 0, 0, 0}, {1, 0, 0, 1}, {1, 0, 1, 0}, {1, 0, 1, 1},
   {1, 1, 0, 0}, {1, 1, 0, 1}, {1, 1, 1, 0}, {1, 1, 1, 1}};

// a has 64 elements. b is a bitvector of 64 dimensions.
float dot(float * restrict a, uint64_t b) {
     int i;
    float sum = 0;
    for(i=0; b && i<64; i+=4,b>>=4) {
        __m128 t0 = _mm_load_ps (a);
        a += 4;
        __m128 t1 = _mm_load_ps (fltmask[b & 15]);
        sum += _mm_cvtss_f32 (_mm_dp_ps (t0, t1, 15));
    }
    return sum;
}

PS。配列は 16 バイトでアラインする必要がaあります。fltmask

PPS。ループでコンパイルするとgcc -std=c99 -msse4 -O2、次のようになります。

.L3:
    movq    %rdx, %rax
    movaps  (%rcx), %xmm1
    shrq    $4, %rdx
    andl    $15, %eax
    addq    $16, %rcx
    addl    $4, %r8d
    salq    $4, %rax
    testq   %rdx, %rdx
    dpps    $15, (%r9,%rax), %xmm1
    addss   %xmm1, %xmm0
    jne .L13

-O3もちろん、展開した状態で。

于 2013-10-11T05:24:12.453 に答える
0

次のようにブランチを削除できます。

for(int i=0; b && i<64; i++, b>>=1)
    sum += a[i] * (b & 1);

これによりマルチが追加されますが、少なくともパイプラインが破壊されることはありません。

分岐を制御するもう 1 つの方法は、分岐をそのまま使用することですが、コンパイラ マクロも使用します。likely(if ...)gccではマクロだと思います。ブランチを使用しますが、そのようにして、ブランチがより頻繁に実行されることをコンパイラに伝え、gcc はより最適化します。

実行できる別の最適化は、内積の「キャッシュ」です。したがって、内積を計算する関数を使用する代わりに、最初に積を 0 に初期化した変数を使用します。また、ベクトルの要素を挿入/削除/更新するたびに、結果を保持する変数も更新します。

于 2013-10-12T16:24:25.587 に答える