c - 在禁用优化的情况下，演示代码无法显示 4 倍快的 SIMD 速度

Question

我试图了解使用 SIMD 矢量化的好处，并编写了一个简单的演示程序代码来查看利用矢量化 (SIMD) 的算法相对于另一个算法的速度增益是多少。以下是 2 种算法:

Alg_A - 不支持 vector :

#include <stdio.h>

#define SIZE 1000000000

int main() {
    printf("Algorithm with NO vector support\n");

    int a[] = {1, 2, 3, 4};
    int b[] = {5, 6, 7, 8};
    int i = 0;

    printf("Running loop %d times\n", SIZE);
    for (; i < SIZE; i++) {
        a[0] = a[0] + b[0];
        a[1] = a[1] + b[1];
        a[2] = a[2] + b[2];
        a[3] = a[3] + b[3];
    }

    printf("A: [%d %d %d %d]\n", a[0], a[1], a[2], a[3]);
}

Alg_B - 支持 vector :

#include <stdio.h>

#define SIZE 1000000000

typedef int v4_i __attribute__ ((vector_size(16)));
union Vec4 {
    v4_i v;
    int i[4];
};

int main() {
    printf("Algorithm with vector support\n\n");

    union Vec4 a, b;
    a.i[0] = 1, a.i[1] = 2, a.i[2] = 3, a.i[3] = 4;
    b.i[0] = 5, b.i[1] = 5, b.i[2] = 7, b.i[3] = 8;
    int i = 0;
    printf("Running loop %d times\n", SIZE);
    for (; i < SIZE; i++) {
        a.v = a.v + b.v;
    }

    printf("A: [%d %d %d %d]\n", a.i[0], a.i[1], a.i[2], a.i[3]);
}

编译完成如下:

算法A:

gcc -ggdb -mno-sse -mno-sse2 -mno-sse3 -mno-sse4 -mno-sse4.1 -mno-sse4.2 -c non_vector_support.c
gcc non_vector_support.o -o non_vector_support

算法B:

gcc -ggdb -c vector_support.c
gcc vector_support.o -o vector_support

查看两种算法的生成代码，我可以看到编译没有做任何像“自动矢量化”这样的技巧(例如，看不到 xmm 寄存器):

算法A:

    for (; i < SIZE; i++) {
  74:   eb 30                   jmp    a6 <main+0xa6>
        a[0] = a[0] + b[0];
  76:   8b 55 d8                mov    -0x28(%rbp),%edx
  79:   8b 45 e8                mov    -0x18(%rbp),%eax
  7c:   01 d0                   add    %edx,%eax
  7e:   89 45 d8                mov    %eax,-0x28(%rbp)
        a[1] = a[1] + b[1];
  81:   8b 55 dc                mov    -0x24(%rbp),%edx
  84:   8b 45 ec                mov    -0x14(%rbp),%eax
  87:   01 d0                   add    %edx,%eax
  89:   89 45 dc                mov    %eax,-0x24(%rbp)
        a[2] = a[2] + b[2];
  8c:   8b 55 e0                mov    -0x20(%rbp),%edx
  8f:   8b 45 f0                mov    -0x10(%rbp),%eax
  92:   01 d0                   add    %edx,%eax
  94:   89 45 e0                mov    %eax,-0x20(%rbp)
        a[3] = a[3] + b[3];
  97:   8b 55 e4                mov    -0x1c(%rbp),%edx
  9a:   8b 45 f4                mov    -0xc(%rbp),%eax
  9d:   01 d0                   add    %edx,%eax
  9f:   89 45 e4                mov    %eax,-0x1c(%rbp)
    int a[] = {1, 2, 3, 4};
    int b[] = {5, 6, 7, 8};
    int i = 0;

    printf("Running loop %d times\n", SIZE);
    for (; i < SIZE; i++) {
  a2:   83 45 d4 01             addl   $0x1,-0x2c(%rbp)
  a6:   81 7d d4 ff c9 9a 3b    cmpl   $0x3b9ac9ff,-0x2c(%rbp)
  ad:   7e c7                   jle    76 <main+0x76>
        a[1] = a[1] + b[1];
        a[2] = a[2] + b[2];
        a[3] = a[3] + b[3];
    }

    printf("A: [%d %d %d %d]\n", a[0], a[1], a[2], a[3]);

算法B:

    for (; i < SIZE; i++) {
  74:   eb 16                   jmp    8c <main+0x8c>
        a.v = a.v + b.v;
  76:   66 0f 6f 4d d0          movdqa -0x30(%rbp),%xmm1
  7b:   66 0f 6f 45 e0          movdqa -0x20(%rbp),%xmm0
  80:   66 0f fe c1             paddd  %xmm1,%xmm0
  84:   0f 29 45 d0             movaps %xmm0,-0x30(%rbp)
    union Vec4 a, b;
    a.i[0] = 1, a.i[1] = 2, a.i[2] = 3, a.i[3] = 4;
    b.i[0] = 5, b.i[1] = 5, b.i[2] = 7, b.i[3] = 8;
    int i = 0;
    printf("Running loop %d times\n", SIZE);
    for (; i < SIZE; i++) {
  88:   83 45 cc 01             addl   $0x1,-0x34(%rbp)
  8c:   81 7d cc ff c9 9a 3b    cmpl   $0x3b9ac9ff,-0x34(%rbp)
  93:   7e e1                   jle    76 <main+0x76>
        a.v = a.v + b.v;
    }

    printf("A: [%d %d %d %d]\n", a.i[0], a.i[1], a.i[2], a.i[3]);

当我运行这些程序时，我希望看到速度提高 4 倍，但是对于这种数据大小，增益似乎平均为 =~ 1s，如果将 SIZE 增加到 8000000000 左右，增益是=~ 5s。这是上面代码中值的时间:

算法A:

Algorithm with NO vector support
Running loop 1000000000 times
A: [705032705 1705032706 -1589934589 -589934588]

real    0m3.279s
user    0m3.282s
sys     0m0.000s

算法B:

支持 vector 的算法

Running loop 1000000000 times
A: [705032705 705032706 -1589934589 -589934588]

real    0m2.609s
user    0m2.607s
sys     0m0.004s

计算与循环相关的开销。我针对给定的 SIZE 运行了一个空循环，并获得了 =~ 2.2s 的平均值。这使我的平均速度提高了大约 2.5 倍。

我是否遗漏了代码或编译行中的某些内容？或者，这是否应该是正确的，在这种情况下，如果我在每次迭代中利用 4 个数据点和 1 条指令，为什么有人会知道为什么速度没有提高 4 倍？

提前致谢。

Answer 1

我在下面整理了一些示例代码，以说明您如何看待 SIMD 与标量代码的优势。示例代码有点做作，但要注意的要点是循环中需要有足够的算术运算以减轻加载/存储延迟和循环开销 - 与初始实验一样，单个添加操作是不够的.

此示例将 32 位 int 数据的吞吐量提高了大约 4 倍。 SIMD 循环有两种版本:一种是没有展开的简单循环，另一种是有 2 次展开的替代循环。正如预期的那样，展开的循环要快一些。

#include <assert.h>
#include <stdio.h>
#include <stdint.h>
#include <stdlib.h>
#include <sys/time.h>   // gettimeofday
#include <smmintrin.h>  // SSE 4.1

static void foo_scalar(uint32_t *a, const uint32_t *b, const uint32_t *c, size_t n)
{
    for (size_t i = 0; i < n; ++i)
    {
        a[i] = (b[i] + c[i] + 1) / 2;
    }
}

static void foo_simd(uint32_t *a, const uint32_t *b, const uint32_t *c, size_t n)
{
    size_t i;

#ifndef UNROLL
    for (i = 0; i <= n - 4; i += 4)
    {
        __m128i vb = _mm_loadu_si128((__m128i *)&b[i]);
        __m128i vc = _mm_loadu_si128((__m128i *)&c[i]);
        __m128i v = _mm_add_epi32(vb, vc);
        v = _mm_add_epi32(v, _mm_set1_epi32(1));
        v = _mm_srli_epi32(v, 1);
        _mm_storeu_si128((__m128i *)&a[i], v);
    }
#else
    for (i = 0; i <= n - 8; i += 8)
    {
        __m128i vb0 = _mm_loadu_si128((__m128i *)&b[i]);
        __m128i vb1 = _mm_loadu_si128((__m128i *)&b[i + 4]);
        __m128i vc0 = _mm_loadu_si128((__m128i *)&c[i]);
        __m128i vc1 = _mm_loadu_si128((__m128i *)&c[i + 4]);
        __m128i v0 = _mm_add_epi32(vb0, vc0);
        __m128i v1 = _mm_add_epi32(vb1, vc1);
        v0 = _mm_add_epi32(v0, _mm_set1_epi32(1));
        v1 = _mm_add_epi32(v1, _mm_set1_epi32(1));
        v0 = _mm_srli_epi32(v0, 1);
        v1 = _mm_srli_epi32(v1, 1);
        _mm_storeu_si128((__m128i *)&a[i], v0);
        _mm_storeu_si128((__m128i *)&a[i + 4], v1);
    }
#endif
    foo_scalar(&a[i], &b[i], &c[i], n - i);
}

int main(int argc, char *argv[])
{
    const size_t kLoops = 100000;
    size_t n = 2 * 1024;
    struct timeval t0, t1;
    double t_scalar_ms, t_simd_ms;

    if (argc > 1)
    {
        n = atoi(argv[1]);
    }

    printf("kLoops = %zu, n = %zu\n", kLoops, n);

    uint32_t * a_scalar = malloc(n * sizeof(uint32_t));
    uint32_t * a_simd = malloc(n * sizeof(uint32_t));
    uint32_t * b = malloc(n * sizeof(uint32_t));
    uint32_t * c = malloc(n * sizeof(uint32_t));

    for (size_t i = 0; i < n; ++i)
    {
        a_scalar[i] = a_simd[i] = 0;
        b[i] = rand();
        c[i] = rand();
    }

    gettimeofday(&t0, NULL);
    for (size_t k = 0; k < kLoops; ++k)
    {
        foo_scalar(a_scalar, b, c, n);
    }
    gettimeofday(&t1, NULL);
    t_scalar_ms = ((double)(t1.tv_sec - t0.tv_sec) + (double)(t1.tv_usec - t0.tv_usec) * 1.0e-6) * 1.0e3;

    gettimeofday(&t0, NULL);
    for (size_t k = 0; k < kLoops; ++k)
    {
        foo_simd(a_simd, b, c, n);
    }
    gettimeofday(&t1, NULL);
    t_simd_ms = ((double)(t1.tv_sec - t0.tv_sec) + (double)(t1.tv_usec - t0.tv_usec) * 1.0e-6) * 1.0e3;

    int64_t sum_scalar = 0, sum_simd = 0;
    for (size_t i = 0; i < n; ++i)
    {
        sum_scalar += a_scalar[i];
        sum_simd += a_simd[i];
    }
    assert(sum_scalar == sum_simd);

    printf("t_scalar = %8g ms = %8g ns / point\n", t_scalar_ms, t_scalar_ms / (kLoops * n) * 1e6);
    printf("t_simd   = %8g ms = %8g ns / point\n", t_simd_ms, t_simd_ms / (kLoops * n) * 1e6);
    printf("Speed-up = %2.1fx\n",  t_scalar_ms / t_simd_ms);

    return 0;
}

编译并运行(没有 SIMD 循环展开):

$ gcc-4.8 -fno-tree-vectorize -std=gnu99 -Wall gros_lalo.c -O3 -msse4.1 && ./a.out
kLoops = 100000, n = 2048
t_scalar =  122.668 ms = 0.598965 ns / point
t_simd   =   33.785 ms = 0.164966 ns / point
Speed-up = 3.6x

编译并运行(2x SIMD 循环展开):

$ gcc-4.8 -fno-tree-vectorize -std=gnu99 -Wall gros_lalo.c -O3 -msse4.1 -DUNROLL && ./a.out
kLoops = 100000, n = 2048
t_scalar =  121.897 ms =   0.5952 ns / point
t_simd   =    29.07 ms = 0.141943 ns / point
Speed-up = 4.2x

看看生成的代码很有意思:

标量:

    xorl    %ecx, %ecx
    .align 4
L10:
    movl    0(%rbp,%rcx,4), %esi
    addl    (%rbx,%rcx,4), %esi
    addl    $1, %esi
    shrl    %esi
    movl    %esi, (%r15,%rcx,4)
    addq    $1, %rcx
    cmpq    %r12, %rcx
    jne L10

SIMD(不展开):

    xorl    %ecx, %ecx
    xorl    %eax, %eax
    .align 4
L18:
    movdqu  0(%rbp,%rcx), %xmm2
    addq    $4, %rax
    movdqu  (%rbx,%rcx), %xmm1
    paddd   %xmm2, %xmm1
    paddd   %xmm3, %xmm1
    psrld   $1, %xmm1
    movdqu  %xmm1, (%r14,%rcx)
    addq    $16, %rcx
    cmpq    %r9, %rax
    jbe L18

SIMD(2 倍展开):

    xorl    %edx, %edx
    xorl    %ecx, %ecx
    .align 4
L18:
    movdqu  0(%rbp,%rdx), %xmm5
    addq    $8, %rcx
    movdqu  (%r11,%rdx), %xmm4
    movdqu  (%rbx,%rdx), %xmm2
    movdqu  (%r10,%rdx), %xmm1
    paddd   %xmm5, %xmm2
    paddd   %xmm4, %xmm1
    paddd   %xmm3, %xmm2
    paddd   %xmm3, %xmm1
    psrld   $1, %xmm2
    psrld   $1, %xmm1
    movdqu  %xmm2, 0(%r13,%rdx)
    movdqu  %xmm1, (%rax,%rdx)
    addq    $32, %rdx
    cmpq    %r15, %rcx
    jbe L18

请注意，前两个循环中的指令数量相似，但 SIMD 循环当然每次迭代处理四个元素，而标量循环每次迭代仅处理一个元素。对于第三个展开的循环，我们有更多的指令，但我们每次迭代处理八个元素 - 请注意，相对于没有循环展开的 SIMD 循环，循环管理指令的比例已经减少。

时序数据是使用 2.6 GHz Core i7 Haswell CPU 在 Mac OS X 10.10 上使用 gcc 4.8 收集的。然而，性能结果在任何当前合理的 x86 CPU 上应该是相似的。