9

GCC (私は 4.8.4 を使用しています) にwhile一番下の関数のループを完全に展開するように指示する方法はありますか?つまり、このループをピールしますか? ループの反復回数は、コンパイル時にわかっています: 58。

最初に私が試したことを説明しましょう。

GAS出力を確認することにより：

gcc -fpic -O2 -S GEPDOT.c

XMM0 ～ XMM11 の 12 本のレジスタを使用します。フラグ-funroll-loopsを gcc に渡すと、次のようになります。

gcc -fpic -O2 -funroll-loops -S GEPDOT.c

ループは 2 回だけ展開されます。GCC 最適化オプションを確認しました。GCC はそれ-funroll-loopsもオンにすると言い-frename-registersます。そのため、GCC がループをアンロールするとき、レジスタ割り当てのための以前の選択は、「残りの」レジスタを使用することです。しかし、XMM12 から XMM15 までは 4 つしか残っていないため、GCC は最高でも 2 回しか展開できません。利用可能な XMM レジスタが 16 個ではなく 48 個ある場合、GCC は問題なく while ループを 4 回アンロールします。

それでも私は別の実験をしました。最初に while ループを手動で 2 回展開し、 function を取得しましたGEPDOT_2。そしたら全然違いがない

gcc -fpic -O2 -S GEPDOT_2.c

と

gcc -fpic -O2 -funroll-loops -S GEPDOT_2.c

GEPDOT_2すでにすべてのレジスタを使用しているため、展開は実行されません。

GCC は、誤った依存関係が導入される可能性を回避するために、レジスタの名前変更を行います。しかし、私の中にはそのような可能性がないことは確かGEPDOTです。あったとしても、それは重要ではありません。ループを自分で展開してみましたが、4 回展開すると 2 回展開するよりも速く、展開しないよりも高速です。もちろん、手動で何度も展開できますが、面倒です。GCCは私のためにこれを行うことができますか? ありがとう。

// C file "GEPDOT.c"
#include <emmintrin.h>

void GEPDOT (double *A, double *B, double *C) {
  __m128d A1_vec = _mm_load_pd(A); A += 2;
  __m128d B_vec = _mm_load1_pd(B); B++;
  __m128d C1_vec = A1_vec * B_vec;
  __m128d A2_vec = _mm_load_pd(A); A += 2;
  __m128d C2_vec = A2_vec * B_vec;
  B_vec = _mm_load1_pd(B); B++;
  __m128d C3_vec = A1_vec * B_vec;
  __m128d C4_vec = A2_vec * B_vec;
  B_vec = _mm_load1_pd(B); B++;
  __m128d C5_vec = A1_vec * B_vec;
  __m128d C6_vec = A2_vec * B_vec;
  B_vec = _mm_load1_pd(B); B++;
  __m128d C7_vec = A1_vec * B_vec;
  A1_vec = _mm_load_pd(A); A += 2;
  __m128d C8_vec = A2_vec * B_vec;
  B_vec = _mm_load1_pd(B); B++;
  int k = 58;
  /* can compiler unroll the loop completely (i.e., peel this loop)? */
  while (k--) {
    C1_vec += A1_vec * B_vec;
    A2_vec = _mm_load_pd(A); A += 2;
    C2_vec += A2_vec * B_vec;
    B_vec = _mm_load1_pd(B); B++;
    C3_vec += A1_vec * B_vec;
    C4_vec += A2_vec * B_vec;
    B_vec = _mm_load1_pd(B); B++;
    C5_vec += A1_vec * B_vec;
    C6_vec += A2_vec * B_vec;
    B_vec = _mm_load1_pd(B); B++;
    C7_vec += A1_vec * B_vec;
    A1_vec = _mm_load_pd(A); A += 2;
    C8_vec += A2_vec * B_vec;
    B_vec = _mm_load1_pd(B); B++;
    }
  C1_vec += A1_vec * B_vec;
  A2_vec = _mm_load_pd(A);
  C2_vec += A2_vec * B_vec;
  B_vec = _mm_load1_pd(B); B++;
  C3_vec += A1_vec * B_vec;
  C4_vec += A2_vec * B_vec;
  B_vec = _mm_load1_pd(B); B++;
  C5_vec += A1_vec * B_vec;
  C6_vec += A2_vec * B_vec;
  B_vec = _mm_load1_pd(B);
  C7_vec += A1_vec * B_vec;
  C8_vec += A2_vec * B_vec;
  /* [write-back] */
  A1_vec = _mm_load_pd(C); C1_vec = A1_vec - C1_vec;
  A2_vec = _mm_load_pd(C + 2); C2_vec = A2_vec - C2_vec;
  A1_vec = _mm_load_pd(C + 4); C3_vec = A1_vec - C3_vec;
  A2_vec = _mm_load_pd(C + 6); C4_vec = A2_vec - C4_vec;
  A1_vec = _mm_load_pd(C + 8); C5_vec = A1_vec - C5_vec;
  A2_vec = _mm_load_pd(C + 10); C6_vec = A2_vec - C6_vec;
  A1_vec = _mm_load_pd(C + 12); C7_vec = A1_vec - C7_vec;
  A2_vec = _mm_load_pd(C + 14); C8_vec = A2_vec - C8_vec;
  _mm_store_pd(C,C1_vec); _mm_store_pd(C + 2,C2_vec);
  _mm_store_pd(C + 4,C3_vec); _mm_store_pd(C + 6,C4_vec);
  _mm_store_pd(C + 8,C5_vec); _mm_store_pd(C + 10,C6_vec);
  _mm_store_pd(C + 12,C7_vec); _mm_store_pd(C + 14,C8_vec);
  }

更新 1

@ user3386109 のコメントのおかげで、この質問を少し広げたいと思います。@ user3386109 は非常に良い質問をします。実際、スケジュールする並列命令が非常に多い場合に、最適なレジスタ割り当てを行うコンパイラの能力に疑問があります。

個人的には、最初にasmインライン アセンブリでループ本体 (HPC の鍵となる) をコーディングし、それを何度でも複製するのが確実な方法だと思います。今年の初めに人気のない投稿がありました: inline assembly . ループの反復回数 j は関数の引数であるため、コンパイル時には不明であるため、コードは少し異なります。その場合、ループを完全に展開することはできないので、アセンブリ コードを 2 回だけ複製し、ループをラベルに変換してジャンプしました。作成したアセンブリの結果として得られるパフォーマンスは、コンパイラが生成したアセンブリよりも約 5% 高いことがわかりました。これは、コンパイラが予想される最適な方法でレジスタを割り当てることができないことを示唆している可能性があります。

私は (そして今も) アセンブリ コーディングの初心者だったので、x86 アセンブリについて少し学ぶための良いケース スタディになりました。GEPDOTしかし、長い目で見れば、私はアセンブリの割合が大きいコードを書く傾向はありません。主に次の3つの理由があります。

  1. asmインライン アセンブリは、移植性がないと批判されています。理由はわかりませんが。おそらく、マシンごとに異なるレジスタが上書きされているためでしょうか?
  2. コンパイラも良くなっています。したがって、コンパイラーが適切な出力を生成するのを支援するために、アルゴリズムの最適化と C コーディングの習慣を改善することを引き続き好みます。
  3. 最後の理由はもっと重要です。反復回数は必ずしも 58 回とは限りません。高性能の行列分解サブルーチンを開発しています。キャッシュ ブロック ファクターのnb場合、反復回数は になりますnb-2。nb以前の投稿で行ったように、関数の引数として指定するつもりはありません。これはマシン固有のパラメータで、マクロとして定義されます。したがって、反復回数はコンパイル時にわかりますが、マシンごとに異なる場合があります。さまざまなnb. したがって、ループを剥がすようにコンパイラに指示するだけの方法があれば、それは素晴らしいことです。

高性能でありながらポータブルなライブラリを作成する経験も共有していただければ幸いです。

4

2 に答える 2

4

これは答えではありませんが、GCC で行列の乗算をベクトル化しようとしている他の人にとっては興味深いかもしれません。

以下では、 cは行優先順の 4×4 行列、 aは列優先順の 4 行n列の行列 (転置)、bは 4 列n行の行列であると仮定します。 -major order で、計算する操作はc = a × b + cです。ここで、× は行列の乗算を表します。

これを達成するための単純な関数は次のとおりです。

void slow_4(double       *c,
            const double *a,
            const double *b,
            size_t        n)
{
    size_t row, col, i;

    for (row = 0; row < 4; row++)
        for (col = 0; col < 4; col++)
            for (i = 0; i < n; i++)
                c[4*row+col] += a[4*i+row] * b[4*i+col];
}

GCC は、以下を使用して SSE2/SSE3 用の非常に優れたコードを生成します。

#if defined(__SSE2__) || defined(__SSE3__)

typedef  double  vec2d  __attribute__((vector_size (2 * sizeof (double))));

void fast_4(vec2d       *c,
            const vec2d *a,
            const vec2d *b,
            size_t       n)
{
    const vec2d *const b_end = b + 2L * n;

    vec2d s00 = c[0];
    vec2d s02 = c[1];
    vec2d s10 = c[2];
    vec2d s12 = c[3];
    vec2d s20 = c[4];
    vec2d s22 = c[5];
    vec2d s30 = c[6];
    vec2d s32 = c[7];

    while (b < b_end) {
        const vec2d b0 = b[0];
        const vec2d b2 = b[1];
        const vec2d a0 = { a[0][0], a[0][0] };
        const vec2d a1 = { a[0][1], a[0][1] };
        const vec2d a2 = { a[1][0], a[1][0] };
        const vec2d a3 = { a[1][1], a[1][1] };
        s00 += a0 * b0;
        s10 += a1 * b0;
        s20 += a2 * b0;
        s30 += a3 * b0;
        s02 += a0 * b2;
        s12 += a1 * b2;
        s22 += a2 * b2;
        s32 += a3 * b2;
        b += 2;
        a += 2;
    }

    c[0] = s00;
    c[1] = s02;
    c[2] = s10;
    c[3] = s12;
    c[4] = s20;
    c[5] = s22;
    c[6] = s30;
    c[7] = s32; 
}

#endif

AVX の場合、GCC は

#if defined(__AVX__) || defined(__AVX2__)

typedef  double  vec4d  __attribute__((vector_size (4 * sizeof (double))));

void fast_4(vec4d       *c,
            const vec4d *a,
            const vec4d *b,
            size_t       n)
{
    const vec4d *const b_end = b + n;

    vec4d s0 = c[0];
    vec4d s1 = c[1];
    vec4d s2 = c[2];
    vec4d s3 = c[3];

    while (b < b_end) {
        const vec4d bc = *(b++);
        const vec4d ac = *(a++);
        const vec4d a0 = { ac[0], ac[0], ac[0], ac[0] };
        const vec4d a1 = { ac[1], ac[1], ac[1], ac[1] };
        const vec4d a2 = { ac[2], ac[2], ac[2], ac[2] };
        const vec4d a3 = { ac[3], ac[3], ac[3], ac[3] };
        s0 += a0 * bc;
        s1 += a1 * bc;
        s2 += a2 * bc;
        s3 += a3 * bc;
    }

    c[0] = s0;
    c[1] = s1;
    c[2] = s2;
    c[3] = s3;
}

#endif

-O2 -march=x86-64 -mtune=generic -msse3gcc-4.8.4 ( )を使用して生成されたアセンブリの SSE3 バージョンは、基本的には

fast_4:
        salq    $5, %rcx
        movapd  (%rdi), %xmm13
        addq    %rdx, %rcx
        cmpq    %rcx, %rdx
        movapd  16(%rdi), %xmm12
        movapd  32(%rdi), %xmm11
        movapd  48(%rdi), %xmm10
        movapd  64(%rdi), %xmm9
        movapd  80(%rdi), %xmm8
        movapd  96(%rdi), %xmm7
        movapd  112(%rdi), %xmm6
        jnb     .L2
.L3:
        movddup (%rsi), %xmm5
        addq    $32, %rdx
        movapd  -32(%rdx), %xmm1
        addq    $32, %rsi
        movddup -24(%rsi), %xmm4
        movapd  %xmm5, %xmm14
        movddup -16(%rsi), %xmm3
        movddup -8(%rsi), %xmm2
        mulpd   %xmm1, %xmm14
        movapd  -16(%rdx), %xmm0
        cmpq    %rdx, %rcx
        mulpd   %xmm0, %xmm5
        addpd   %xmm14, %xmm13
        movapd  %xmm4, %xmm14
        mulpd   %xmm0, %xmm4
        addpd   %xmm5, %xmm12
        mulpd   %xmm1, %xmm14
        addpd   %xmm4, %xmm10
        addpd   %xmm14, %xmm11
        movapd  %xmm3, %xmm14
        mulpd   %xmm0, %xmm3
        mulpd   %xmm1, %xmm14
        mulpd   %xmm2, %xmm0
        addpd   %xmm3, %xmm8
        mulpd   %xmm2, %xmm1
        addpd   %xmm14, %xmm9
        addpd   %xmm0, %xmm6
        addpd   %xmm1, %xmm7
        ja      .L3
.L2:
        movapd  %xmm13, (%rdi)
        movapd  %xmm12, 16(%rdi)
        movapd  %xmm11, 32(%rdi)
        movapd  %xmm10, 48(%rdi)
        movapd  %xmm9, 64(%rdi)
        movapd  %xmm8, 80(%rdi)
        movapd  %xmm7, 96(%rdi)
        movapd  %xmm6, 112(%rdi)
        ret

生成されたアセンブリ ( -O2 -march=x86-64 -mtune=generic -mavx) の AVX バージョンは基本的に

fast_4:
        salq       $5, %rcx
        vmovapd    (%rdi), %ymm5
        addq       %rdx, %rcx
        vmovapd    32(%rdi), %ymm4
        cmpq       %rcx, %rdx
        vmovapd    64(%rdi), %ymm3
        vmovapd    96(%rdi), %ymm2
        jnb        .L2
.L3:
        addq       $32, %rsi
        vmovapd    -32(%rsi), %ymm1
        addq       $32, %rdx
        vmovapd    -32(%rdx), %ymm0
        cmpq       %rdx, %rcx
        vpermilpd  $0, %ymm1, %ymm6
        vperm2f128 $0, %ymm6, %ymm6, %ymm6
        vmulpd     %ymm0, %ymm6, %ymm6
        vaddpd     %ymm6, %ymm5, %ymm5
        vpermilpd  $15, %ymm1, %ymm6
        vperm2f128 $0, %ymm6, %ymm6, %ymm6
        vmulpd     %ymm0, %ymm6, %ymm6
        vaddpd     %ymm6, %ymm4, %ymm4
        vpermilpd  $0, %ymm1, %ymm6
        vpermilpd  $15, %ymm1, %ymm1
        vperm2f128 $17, %ymm6, %ymm6, %ymm6
        vperm2f128 $17, %ymm1, %ymm1, %ymm1
        vmulpd     %ymm0, %ymm6, %ymm6
        vmulpd     %ymm0, %ymm1, %ymm0
        vaddpd     %ymm6, %ymm3, %ymm3
        vaddpd     %ymm0, %ymm2, %ymm2
        ja         .L3
.L2:
        vmovapd    %ymm5, (%rdi)
        vmovapd    %ymm4, 32(%rdi)
        vmovapd    %ymm3, 64(%rdi)
        vmovapd    %ymm2, 96(%rdi)
        vzeroupper
        ret

レジスターのスケジューリングは最適ではないと思いますが、それほどひどいものでもありません。この時点で手動で最適化しようとせずに、個人的には上記に満足しています。

Core i5-4200U プロセッサ (AVX2 対応) では、上記の関数の高速バージョンは、SSE3 では 1843 CPU サイクル (中央値)、AVX2 では 1248 サイクルで 2 つの 4×256 行列の積を計算します。これは、行列エントリごとに 1.8 および 1.22 サイクルになります。比較のため、ベクトル化されていない低速バージョンでは、行列のエントリごとに約 11 サイクルかかります。

(サイクル カウントは中央値です。つまり、テストの半分の方が高速でした。約 10 万回程度の繰り返しで大まかなベンチマークを実行しただけなので、これらの数値は割り引いて考えてください。)

この CPU では、キャッシュ効果により、4×512 マトリックス サイズの AVX2 はエントリあたり 1.2 サイクルのままですが、4×1024 では 1.4 に、4×4096 では 1.5 に、4×8192 では 1.8 に、そして4×65536 でエントリあたり 2.2 サイクル。SSE3 バージョンは、4×3072 までエントリあたり 1.8 サイクルのままであり、その時点で速度が低下し始めます。4×65536 では、1 エントリあたり約 2.2 サイクルです。この (ラップトップ!) CPU は、現時点ではキャッシュ帯域幅が制限されていると思います。

于 2016-03-20T13:09:40.610 に答える
3

オプティマイザーのパラメーターを微調整してみてください。

gcc -O3 -funroll-loops --param max-completely-peeled-insns=1000 --param max-completely-peel-times=100

これでうまくいくはずです。

于 2016-03-20T09:40:18.393 に答える