37

我正在研究矢量化对程序性能的影响。在这方面，我编写了以下代码：

#include <stdio.h>
#include <sys/time.h>
#include <stdlib.h>

#define LEN 10000000

int main(){

    struct timeval stTime, endTime;

    double* a = (double*)malloc(LEN*sizeof(*a));
    double* b = (double*)malloc(LEN*sizeof(*b));
    double* c = (double*)malloc(LEN*sizeof(*c));

    int k;
    for(k = 0; k < LEN; k++){
        a[k] = rand();
        b[k] = rand();
    }

    gettimeofday(&stTime, NULL);

    for(k = 0; k < LEN; k++)
        c[k] = a[k] * b[k];

    gettimeofday(&endTime, NULL);

    FILE* fh = fopen("dump", "w");
    for(k = 0; k < LEN; k++)
        fprintf(fh, "c[%d] = %f\t", k, c[k]);
    fclose(fh);

    double timeE = (double)(endTime.tv_usec + endTime.tv_sec*1000000 - stTime.tv_usec - stTime.tv_sec*1000000);

    printf("Time elapsed: %f\n", timeE);

    return 0;
}

在这段代码中，我只是对两个向量进行初始化和相乘。结果保存在 vector 中c。我主要感兴趣的是矢量化以下循环的效果：

for(k = 0; k < LEN; k++)
    c[k] = a[k] * b[k];

我使用以下两个命令编译代码：

1) icc -O2 TestSMID.c -o TestSMID -no-vec -no-simd
2) icc -O2 TestSMID.c -o TestSMID -vec-report2

我希望看到性能改进，因为第二个命令成功地矢量化了循环。但是，我的研究表明，当循环被矢量化时，性能并没有提高。

我可能在这里错过了一些东西，因为我对这个话题不是很熟悉。所以，如果我的代码有问题，请告诉我。

在此先感谢您的帮助。

PS：我使用的是 Mac OSX，所以不需要对齐数据，因为所有分配的内存都是 16 字节对齐的。

编辑：首先感谢大家的评论和回答。我想到了@Mysticial 提出的答案，这里还有一些应该提到的点。首先，正如@Vinska 提到的，c[k]=a[k]*b[k]不只需要一个周期。除了循环索引增量和进行比较以确保k小于 之外LEN，还需要执行其他操作来执行操作。看一下编译器生成的汇编代码，可以看出一个简单的乘法需要不止一个周期。矢量化版本如下所示：

L_B1.9:                         # Preds L_B1.8
        movq      %r13, %rax                                    #25.5
        andq      $15, %rax                                     #25.5
        testl     %eax, %eax                                    #25.5
        je        L_B1.12       # Prob 50%                      #25.5
                                # LOE rbx r12 r13 r14 r15 eax
L_B1.10:                        # Preds L_B1.9
        testb     $7, %al                                       #25.5
        jne       L_B1.32       # Prob 10%                      #25.5
                                # LOE rbx r12 r13 r14 r15
L_B1.11:                        # Preds L_B1.10
        movsd     (%r14), %xmm0                                 #26.16
        movl      $1, %eax                                      #25.5
        mulsd     (%r15), %xmm0                                 #26.23
        movsd     %xmm0, (%r13)                                 #26.9
                                # LOE rbx r12 r13 r14 r15 eax
L_B1.12:                        # Preds L_B1.11 L_B1.9
        movl      %eax, %edx                                    #25.5
        movl      %eax, %eax                                    #26.23
        negl      %edx                                          #25.5
        andl      $1, %edx                                      #25.5
        negl      %edx                                          #25.5
        addl      $10000000, %edx                               #25.5
        lea       (%r15,%rax,8), %rcx                           #26.23
        testq     $15, %rcx                                     #25.5
        je        L_B1.16       # Prob 60%                      #25.5
                                # LOE rdx rbx r12 r13 r14 r15 eax
L_B1.13:                        # Preds L_B1.12
        movl      %eax, %eax                                    #25.5
                                # LOE rax rdx rbx r12 r13 r14 r15
L_B1.14:                        # Preds L_B1.14 L_B1.13
        movups    (%r15,%rax,8), %xmm0                          #26.23
        movsd     (%r14,%rax,8), %xmm1                          #26.16
        movhpd    8(%r14,%rax,8), %xmm1                         #26.16
        mulpd     %xmm0, %xmm1                                  #26.23
        movntpd   %xmm1, (%r13,%rax,8)                          #26.9
        addq      $2, %rax                                      #25.5
        cmpq      %rdx, %rax                                    #25.5
        jb        L_B1.14       # Prob 99%                      #25.5
        jmp       L_B1.20       # Prob 100%                     #25.5
                                # LOE rax rdx rbx r12 r13 r14 r15
L_B1.16:                        # Preds L_B1.12
        movl      %eax, %eax                                    #25.5
                                # LOE rax rdx rbx r12 r13 r14 r15
L_B1.17:                        # Preds L_B1.17 L_B1.16
        movsd     (%r14,%rax,8), %xmm0                          #26.16
        movhpd    8(%r14,%rax,8), %xmm0                         #26.16
        mulpd     (%r15,%rax,8), %xmm0                          #26.23
        movntpd   %xmm0, (%r13,%rax,8)                          #26.9
        addq      $2, %rax                                      #25.5
        cmpq      %rdx, %rax                                    #25.5
        jb        L_B1.17       # Prob 99%                      #25.5
                                # LOE rax rdx rbx r12 r13 r14 r15
L_B1.18:                        # Preds L_B1.17
        mfence                                                  #25.5
                                # LOE rdx rbx r12 r13 r14 r15
L_B1.19:                        # Preds L_B1.18
        mfence                                                  #25.5
                                # LOE rdx rbx r12 r13 r14 r15
L_B1.20:                        # Preds L_B1.14 L_B1.19 L_B1.32
        cmpq      $10000000, %rdx                               #25.5
        jae       L_B1.24       # Prob 0%                       #25.5
                                # LOE rdx rbx r12 r13 r14 r15
L_B1.22:                        # Preds L_B1.20 L_B1.22
        movsd     (%r14,%rdx,8), %xmm0                          #26.16
        mulsd     (%r15,%rdx,8), %xmm0                          #26.23
        movsd     %xmm0, (%r13,%rdx,8)                          #26.9
        incq      %rdx                                          #25.5
        cmpq      $10000000, %rdx                               #25.5
        jb        L_B1.22       # Prob 99%                      #25.5
                                # LOE rdx rbx r12 r13 r14 r15
L_B1.24:                        # Preds L_B1.22 L_B1.20

非矢量化版本是：

L_B1.9:                         # Preds L_B1.8
        xorl      %eax, %eax                                    #25.5
                                # LOE rbx r12 r13 r14 r15 eax
L_B1.10:                        # Preds L_B1.10 L_B1.9
        lea       (%rax,%rax), %edx                             #26.9
        incl      %eax                                          #25.5
        cmpl      $5000000, %eax                                #25.5
        movsd     (%r15,%rdx,8), %xmm0                          #26.16
        movsd     8(%r15,%rdx,8), %xmm1                         #26.16
        mulsd     (%r13,%rdx,8), %xmm0                          #26.23
        mulsd     8(%r13,%rdx,8), %xmm1                         #26.23
        movsd     %xmm0, (%rbx,%rdx,8)                          #26.9
        movsd     %xmm1, 8(%rbx,%rdx,8)                         #26.9
        jb        L_B1.10       # Prob 99%                      #25.5
                                # LOE rbx r12 r13 r14 r15 eax

除此之外，处理器不仅仅加载 24 个字节。在每次访问内存时，都会加载一个完整的行（64 字节）。a更重要的是，由于, b, 和所需的内存c是连续的，预取器肯定会有很大帮助并提前加载下一个块。说了这么多，我觉得@Mysticial计算出来的内存带宽太悲观了。

此外，英特尔矢量化指南中提到了使用 SIMD 来提高程序性能的非常简单的加法。因此，看起来我们应该能够为这个非常简单的循环获得一些性能改进。

Edit2：再次感谢您的评论。另外，感谢@Mysticial 示例代码，我终于看到了SIMD 对性能提升的影响。正如 Mysticial 所提到的，问题在于内存带宽。通过为 、 、 选择较小的尺寸a，b并且c适合 L1 缓存，可以看出 SIMD 可以帮助显着提高性能。这是我得到的结果：

icc -O2 -o TestSMIDNoVec -no-vec TestSMID2.c: 17.34 sec

icc -O2 -o TestSMIDVecNoUnroll -vec-report2 TestSMID2.c: 9.33 sec

展开循环可以进一步提高性能：

icc -O2 -o TestSMIDVecUnroll -vec-report2 TestSMID2.c -unroll=8: 8.6sec

另外，我应该提一下，当使用-O2.

PS：我的电脑是Macbook Pro core i5 @2.5GHz（双核）

4

4 回答 4

78

这个原始答案早在 2013 年就有效。截至 2017 年硬件，情况已经发生了很大变化，以至于问题和答案都已过时。

请参阅此答案的结尾以了解 2017 年更新。

原始答案（2013）：

因为您受到内存带宽的限制。

虽然矢量化和其他微优化可以提高计算速度，但它们不能提高内存速度。

在您的示例中：

for(k = 0; k < LEN; k++)
    c[k] = a[k] * b[k];

您正在对所有内存进行一次遍历，所做的工作很少。这是最大化您的内存带宽。

因此，无论它如何优化（矢量化、展开等），它都不会变得更快。

2013 年的典型台式机具有大约10 GB/s的内存带宽*。
您的循环触及24 bytes/iteration。

如果没有矢量化，现代 x64 处理器可能每个周期可以进行大约 1 次迭代*。

假设您以 4 GHz 运行：

  • (4 * 10^9) * 24 bytes/iteration = 96 GB/s

这几乎是内存带宽的 10 倍——没有矢量化。

*毫不奇怪，因为我没有给出引用，所以一些人怀疑我上面给出的数字。好吧，从经验来看，这些都是我的想法。所以这里有一些基准来证明它。

循环迭代可以以 1 个循环/迭代的速度运行：

LEN如果我们减少它以适应缓存，我们可以摆脱内存瓶颈。
（我在 C++ 中对此进行了测试，因为它更容易。但它没有区别。）

#include <iostream>
#include <time.h>
using std::cout;
using std::endl;

int main(){
    const int LEN = 256;

    double *a = (double*)malloc(LEN*sizeof(*a));
    double *b = (double*)malloc(LEN*sizeof(*a));
    double *c = (double*)malloc(LEN*sizeof(*a));

    int k;
    for(k = 0; k < LEN; k++){
        a[k] = rand();
        b[k] = rand();
    }

    clock_t time0 = clock();

    for (int i = 0; i < 100000000; i++){
        for(k = 0; k < LEN; k++)
            c[k] = a[k] * b[k];
    }

    clock_t time1 = clock();
    cout << (double)(time1 - time0) / CLOCKS_PER_SEC << endl;
}
  • 处理器：英特尔酷睿 i7 2600K @ 4.2 GHz
  • 编译器：Visual Studio 2012
  • 时间：6.55秒

在这个测试中，我仅在6.55秒内运行了 25,600,000,000 次迭代。

  • 6.55 * 4.2 GHz= 27,510,000,000 次循环
  • 27,510,000,000 / 25,600,000,000= 1.074 次循环/迭代

现在，如果您想知道如何做到这一点：

  • 2 个负载
  • 1 家商店
  • 1 乘
  • 递增计数器
  • 比较+分支

全部在一个周期内...

这是因为现代处理器和编译器非常棒。

虽然这些操作中的每一个都有延迟（尤其是乘法），但处理器能够同时执行多个迭代。我的测试机器是 Sandy Bridge 处理器，它能够承受 2x128b 负载、1x128b 存储和 1x256b 矢量 FP 乘以每个周期。如果负载是微融合微指令的内存源操作数，则可能还有另外一两个向量或整数运算。（仅当使用 256b AVX 加载/存储时，2 加载 + 1 存储吞吐量，否则每个周期只有两个总内存操作（最多一个存储））。

查看程序集（为简洁起见，我将省略），编译器似乎展开了循环，从而减少了循环开销。但它并没有完全矢量化它。

内存带宽约为 10 GB/s：

最简单的测试方法是通过memset()：

#include <iostream>
#include <time.h>
using std::cout;
using std::endl;

int main(){
    const int LEN = 1 << 30;    //  1GB

    char *a = (char*)calloc(LEN,1);

    clock_t time0 = clock();

    for (int i = 0; i < 100; i++){
        memset(a,0xff,LEN);
    }

    clock_t time1 = clock();
    cout << (double)(time1 - time0) / CLOCKS_PER_SEC << endl;
}
  • 处理器：英特尔酷睿 i7 2600K @ 4.2 GHz
  • 编译器：Visual Studio 2012
  • 时间：5.811 秒

所以我的机器需要5.811秒才能写入 100 GB 的内存。这大约是17.2 GB/s。

我的处理器处于高端。Nehalem 和 Core 2 代处理器的内存带宽较少。

2017 年 3 月更新：

到 2017 年，事情变得更加复杂。

得益于 DDR4 和四通道内存，单线程不再可能使内存带宽饱和。但带宽问题并不一定会消失。尽管带宽增加了，但处理器内核也得到了改进——而且数量更多。

从数学上讲：

  • 每个核心都有带宽限制X。
  • 主内存的带宽限制为Y.
  • 在旧系统上，X > Y.
  • 在当前的高端系统上，X < Y. 但是X * (# of cores) > Y。

早在 2013 年：Sandy Bridge @ 4 GHz + 双通道 DDR3 @ 1333 MHz

  • 无向量化（8 字节加载/存储）：X = 32 GB/s和Y = ~17 GB/s
  • 矢量化 SSE*（16 字节加载/存储）：X = 64 GB/s和Y = ~17 GB/s

现在在 2017 年：Haswell-E @ 4 GHz + 四通道 DDR4 @ 2400 MHz

  • 无向量化（8 字节加载/存储）：X = 32 GB/s和Y = ~70 GB/s
  • 矢量化 AVX*（32 字节加载/存储）：X = 64 GB/s和Y = ~70 GB/s

（对于 Sandy Bridge 和 Haswell，缓存中的架构限制将带宽限制为大约 16 字节/周期，无论 SIMD 宽度如何。）

因此，如今，单个线程并不总是能够使内存带宽饱和。您需要进行矢量化以达到X. 但是您仍然会达到Y2 个或更多线程的主内存带宽限制。

但是有一件事没有改变，而且很可能很长时间都不会改变：如果总内存带宽不饱和，您将无法在所有内核上运行占用带宽的循环。

于 2013-08-10T07:02:59.160 回答
3

正如 Mysticial 已经描述的，主存带宽限制是这里大缓冲区的瓶颈。解决此问题的方法是重新设计您的处理以在适合缓存的块中工作。（与其将整个 200MiB 的双精度相乘，不如仅乘以 128kiB，然后对其进行处理。因此，使用乘法输出的代码会发现它仍在 L2 缓存中。L2 通常为 256kiB，并且是每个 CPU 内核私有的，关于最近的英特尔设计。）

这种技术称为缓存阻塞或循环平铺。 对于某些算法来说这可能很棘手，但回报是 L2 缓存带宽与主内存带宽之间的差异。

如果您这样做，请确保编译器不会仍在生成流式存储 ( movnt...)。这些写入绕过缓存以避免不适合的数据污染它。该数据的下一次读取将需要触及主存储器。

于 2015-07-05T03:12:45.760 回答
2

编辑：修改了很多答案。另外，请忽略我之前写的关于 Mystical 的答案并不完全正确的大部分内容。尽管如此，我仍然不同意它受到内存的限制，因为尽管进行了各种各样的测试，但我看不到任何原始代码受内存速度限制的迹象。同时，它不断显示出受 CPU 限制的明显迹象。

可能有很多原因。而且由于原因可能非常依赖于硬件，因此我决定不应该根据猜测进行推测。只是概述一下我在以后的测试中遇到的这些事情，我使用了一种更准确可靠的 CPU 时间测量方法和循环 1000 次。我相信这些信息可能会有所帮助。但是请谨慎对待它，因为它取决于硬件。

  • 当使用来自 SSE 系列的指令时，我得到的向量化代码比非向量化代码快 10% 以上。
  • 使用 SSE 系列的矢量化代码和使用 AVX 的矢量化代码或多或少地以相同的性能运行。
  • 使用 AVX 指令时，非矢量化代码运行速度最快——比我尝试过的所有其他方法快 25% 或更多。
  • 在所有情况下，结果都与 CPU 时钟成线性比例。
  • 结果几乎不受内存时钟的影响。
  • 结果受到内存延迟的影响很大 - 远远超过内存时钟，但几乎没有 CPU 时钟对结果的影响那么大。

WRT Mystical 的每个时钟运行近 1 次迭代的示例 - 我没想到 CPU 调度程序会如此高效，并假设每 1.5-2 个时钟周期进行 1 次迭代。但令我惊讶的是，事实并非如此。我确实错了，对此感到抱歉。我自己的 CPU 运行效率更高 - 1.048 个周期/迭代。因此，我可以证明 Mystical 的这部分答案绝对正确。

于 2013-08-10T12:32:58.483 回答
0

以防万一 a[] b[] 和 c[] 争夺 L2 缓存 ::

#include <string.h> /* for memcpy */

 ...

 gettimeofday(&stTime, NULL);

    for(k = 0; k < LEN; k += 4) {
        double a4[4], b4[4], c4[4];
        memcpy(a4,a+k, sizeof a4);
        memcpy(b4,b+k, sizeof b4);
        c4[0] = a4[0] * b4[0];
        c4[1] = a4[1] * b4[1];
        c4[2] = a4[2] * b4[2];
        c4[3] = a4[3] * b4[3];
        memcpy(c+k,c4, sizeof c4);
        }

    gettimeofday(&endTime, NULL);

将运行时间从 98429.000000 减少到 67213.000000；展开循环 8 倍将其减少到 57157.000000。

于 2013-08-10T13:56:25.743 回答