SSE 加法比 + 运算符慢答案

【问题标题】：SSE addition is slower than + operatorSSE 加法比 + 运算符慢
【发布时间】：2018-11-15 05:01:45
【问题描述】：

我试图测试添加 SSE 的速度有多快，但有些地方不对劲。我在堆栈中创建了两个输入数组和一个输出数组，并以两种方式对它们执行加法。它比常规 + 运算符慢。我在这里做错了什么：

#include <iostream>
#include <nmmintrin.h>
#include <chrono>

using namespace std;

#define USE_SSE

typedef chrono::steady_clock::time_point TimeStamp;
typedef chrono::steady_clock Clock;
int main()
{
    const int MAX = 100000 * 4;
    float in1[MAX];
    float in2[MAX];
    float out[MAX];

    memset(out,0,sizeof(float) * MAX);

    for(int i = 0 ; i < MAX ; ++i)
    {
        in1[i] = 1.0f;
        in2[i] = 1.0f;
    }

    TimeStamp start,end;
    start = Clock::now();

    for(int i = 0 ; i < MAX ; i+=4)
    {
#ifdef USE_SSE

        __m128 a = _mm_load_ps(&in1[i]);
        __m128 b = _mm_load_ps(&in2[i]);
        __m128 result = _mm_add_ps(a,b);
        _mm_store_ps(&out[i],result);
#else
        out[0] = in1[0] + in2[0];
        out[1] = in1[1] + in2[1];
        out[2] = in1[2] + in2[2];
        out[3] = in1[3] + in2[3];
#endif
    }


    end = Clock::now();
    double dt = chrono::duration_cast<chrono::nanoseconds>(end-start).count();
    cout<<dt<<endl;

    return 0;
}

这里是内存对齐问题吗？

【问题讨论】：

让 SSE 版本遍历数组而不是标量版本看起来不公平..
Any 体面的优化编译器无论如何都会将非 sse 路径转换为实际使用 SSE。事实上，它不仅会 vectorise 你的循环，它也很可能会展开它。在那种情况下，天真的循环很可能会更快。编译器确实非常聪明，很难（但并非不可能）击败它们。真正确定发生了什么的唯一方法是检查生成的代码 (asm)。还记得在计时代码时打开优化 - 计时非优化代码通常是没有意义的。

标签： c++ x86 sse simd

【解决方案1】：

您的代码中有错误，非 SSE 部分应为：

    out[i+0] = in1[i+0] + in2[i+0];
    out[i+1] = in1[i+1] + in2[i+1];
    out[i+2] = in1[i+2] + in2[i+2];
    out[i+3] = in1[i+3] + in2[i+3];

您应该考虑让您的基准测试运行时间稍长一些，因为测量短时间段是不可靠的。也许，你需要做一些事情来阻止编译器优化你的代码（比如标记out volatile）。始终检查汇编代码，以确保您测量的是什么。

【讨论】：

【解决方案2】：

这是您的基准测试的一个稍微改进的版本，修复了错误，改进了时序，并为标量代码禁用了编译器矢量化（至少对于 gcc 和 clang）：

#include <iostream>
#include <xmmintrin.h>
#include <chrono>

using namespace std;

typedef chrono::steady_clock::time_point TimeStamp;
typedef chrono::steady_clock Clock;

typedef void (*add_func)(const float *in1, const float *in2, volatile float *out, const size_t n);

#ifndef __clang__
__attribute__((optimize("no-tree-vectorize")))
#endif
static void add_scalar(const float *in1, const float *in2, volatile float *out, const size_t n)
{
#ifdef __clang__
    #pragma clang loop vectorize(disable)
#endif
    for (size_t i = 0 ; i < n ; i += 4)
    {
        out[i + 0] = in1[i + 0] + in2[i + 0];
        out[i + 1] = in1[i + 1] + in2[i + 1];
        out[i + 2] = in1[i + 2] + in2[i + 2];
        out[i + 3] = in1[i + 3] + in2[i + 3];
    }
}

static void add_SIMD(const float *in1, const float *in2, volatile float *out, const size_t n)
{
    for (size_t i = 0 ; i < n ; i += 4)
    {
        __m128 a = _mm_loadu_ps(&in1[i]);
        __m128 b = _mm_loadu_ps(&in2[i]);
        __m128 result = _mm_add_ps(a, b);
        _mm_storeu_ps((float *)&out[i], result);
    }
}

static double time_func(const float *in1, const float *in2, volatile float *out, const size_t n, add_func f)
{
    const size_t kLoops = 10000;

    TimeStamp start,end;
    start = Clock::now();

    for (size_t k = 0; k < kLoops; ++k)
    {
        f(in1, in2, out, n);
    }

    end = Clock::now();

    return chrono::duration_cast<chrono::nanoseconds>(end - start).count() / ((double)kLoops * (double)n);
}

int main()
{
    const size_t n = 100000 * 4;
    float *in1 = new float[n];
    float *in2 = new float[n];
    volatile float *out = new float[n]();

    for (size_t i = 0; i < n; ++i)
    {
        in1[i] = (float)i;
        in2[i] = 1.0f;
    }

    double t_scalar = time_func(in1, in2, out, n, add_scalar);
    double t_SIMD = time_func(in1, in2, out, n, add_SIMD);

    cout << "t_scalar = " << t_scalar << " ns / point" << endl;
    cout << "t_SIMD   = " << t_SIMD << " ns / point" << endl;
    cout << "speed-up = " << t_scalar / t_SIMD << "x" << endl;

    delete [] in1;
    delete [] in2;
    delete [] out;

    return 0;
}

Haswell CPU 上的 SSE 提高了大约 1.5 到 1.6 倍。这显然低于可能的 4 倍理论改进，但测试很可能会受到带宽限制，因为您每次迭代仅执行 1 次算术运算，但 2 次加载和 1 次存储：

t_scalar = 0.529723 ns / point
t_SIMD   = 0.329758 ns / point
speed-up = 1.6064x

【讨论】：

我认为 Haswell 有足够的内存带宽用于单核 AVX。部分问题可能是延迟，而不是带宽。
@MSalters：很难说，但这肯定与内存有关 - 如果你将 n 降低到 100 * 4 以便测试完全超出 L1，那么你将获得预期的 4 倍改进。
@MSalters：每个核心时钟周期加载 64 字节 + 存储 32 字节是一个巨大的带宽量。 L1d 可以在一瞬间做到这一点，但显然不能完全维持它。甚至 L2 也是 a[i] = b[i]+c[i] 的一个重要瓶颈。另一个问题是单核不能使 DRAM 饱和，尤其是在具有许多 DRAM 通道的多核 Xeon 中，其中单核带宽更受max_concurrency * latency 的限制，具有更高的离核延迟。 Why is Skylake so much better than Broadwell-E for single-threaded memory throughput?

【解决方案3】：

有时尝试通过添加循环来测试时序来“优化”C++ 代码通常是非常愚蠢的，这就是其中一种情况 :(

你的代码LITERALLY归结为：

int main()
{
    TimeStamp start = Clock::now();
    TimeStamp end = Clock::now();

    double dt = chrono::duration_cast<chrono::nanoseconds>(end-start).count();
    cout<<dt<<endl;

    return 0;
}

编译器并不愚蠢，因此它决定删除您的内部循环（因为输出未使用，因此循环是多余的）。

即使编译器决定保留循环，每次添加都会发出 3 条内存指令。如果您的内存是 1600Mhz，而您的 CPU 是 3200Mhz，那么您的测试只是向您证明您的内存带宽有限。像这样的分析循环没有用，在分析器中测试真实世界的情况总是会更好......

无论如何，回到有问题的循环。让我们将代码放入编译器资源管理器并尝试一些选项......

https://godbolt.org/z/5SJQHb

F0：只是一个基本的、无聊的 C 循环。

for(int i = 0 ; i < MAX ; i++)
{
    out[i] = in1[i] + in2[i];
}

编译器输出这个内部循环：

vmovups ymm0,YMMWORD PTR [rsi+r8*4]
vmovups ymm1,YMMWORD PTR [rsi+r8*4+0x20]
vmovups ymm2,YMMWORD PTR [rsi+r8*4+0x40]
vmovups ymm3,YMMWORD PTR [rsi+r8*4+0x60]
vaddps ymm0,ymm0,YMMWORD PTR [rdx+r8*4]
vaddps ymm1,ymm1,YMMWORD PTR [rdx+r8*4+0x20]
vaddps ymm2,ymm2,YMMWORD PTR [rdx+r8*4+0x40]
vaddps ymm3,ymm3,YMMWORD PTR [rdx+r8*4+0x60]
vmovups YMMWORD PTR [rdi+r8*4],ymm0
vmovups YMMWORD PTR [rdi+r8*4+0x20],ymm1
vmovups YMMWORD PTR [rdi+r8*4+0x40],ymm2
vmovups YMMWORD PTR [rdi+r8*4+0x60],ymm3

展开，每次迭代处理 32xfloats（在 AVX2 中）[+在迭代结束时处理多达 31 个元素的额外代码]

F1：上面的 SSE“优化”循环。 （显然这段代码在循环结束时最多不能处理 3 个元素）

for(int i = 0 ; i < MAX ; i+=4)
{
    __m128 a = _mm_load_ps(&in1[i]);
    __m128 b = _mm_load_ps(&in2[i]);
    __m128 result = _mm_add_ps(a,b);
    _mm_store_ps(&out[i],result);
}

这个输出：

vmovaps xmm0,XMMWORD PTR [rsi+rcx*4]
vaddps xmm0,xmm0,XMMWORD PTR [rdx+rcx*4]
vmovaps XMMWORD PTR [rdi+rcx*4],xmm0
vmovaps xmm0,XMMWORD PTR [rsi+rcx*4+0x10]
vaddps xmm0,xmm0,XMMWORD PTR [rdx+rcx*4+0x10]
vmovaps XMMWORD PTR [rdi+rcx*4+0x10],xmm0
vmovaps xmm0,XMMWORD PTR [rsi+rcx*4+0x20]
vaddps xmm0,xmm0,XMMWORD PTR [rdx+rcx*4+0x20]
vmovaps XMMWORD PTR [rdi+rcx*4+0x20],xmm0
vmovaps xmm0,XMMWORD PTR [rsi+rcx*4+0x30]
vaddps xmm0,xmm0,XMMWORD PTR [rdx+rcx*4+0x30]
vmovaps XMMWORD PTR [rdi+rcx*4+0x30],xmm0

所以编译器已经展开循环，但它已经退回到 SSE（根据要求），所以现在是原始循环性能的一半（不太正确 - 内存带宽将是这里的限制因素）。

F2：您手动展开的 C++ 循环（索引已更正，但仍然无法处理最后 3 个元素）

for(int i = 0 ; i < MAX ; i += 4)
{
    out[i + 0] = in1[i + 0] + in2[i + 0];
    out[i + 1] = in1[i + 1] + in2[i + 1];
    out[i + 2] = in1[i + 2] + in2[i + 2];
    out[i + 3] = in1[i + 3] + in2[i + 3];
}

还有输出：

vmovss xmm0,DWORD PTR [rsi+rax*4]
vaddss xmm0,xmm0,DWORD PTR [rdx+rax*4]
vmovss DWORD PTR [rdi+rax*4],xmm0
vmovss xmm0,DWORD PTR [rsi+rax*4+0x4]
vaddss xmm0,xmm0,DWORD PTR [rdx+rax*4+0x4]
vmovss DWORD PTR [rdi+rax*4+0x4],xmm0
vmovss xmm0,DWORD PTR [rsi+rax*4+0x8]
vaddss xmm0,xmm0,DWORD PTR [rdx+rax*4+0x8]
vmovss DWORD PTR [rdi+rax*4+0x8],xmm0
vmovss xmm0,DWORD PTR [rsi+rax*4+0xc]
vaddss xmm0,xmm0,DWORD PTR [rdx+rax*4+0xc]
vmovss DWORD PTR [rdi+rax*4+0xc],xmm0

好吧，这完全无法矢量化！它一次只处理 1 个添加。好吧，这通常归结为指针别名，所以我将从这里更改函数原型：

void func(float* out, const float* in1, const float* in2, int MAX);

对此：(F4)

void func(
    float* __restrict out, 
    const float* __restrict in1, 
    const float* __restrict in2, 
    int MAX);

现在编译器将输出向量化的内容：

vmovups xmm0,XMMWORD PTR [rsi+rcx*4]
vaddps xmm0,xmm0,XMMWORD PTR [rdx+rcx*4]
vmovups xmm1,XMMWORD PTR [rsi+rcx*4+0x10]
vaddps xmm1,xmm1,XMMWORD PTR [rdx+rcx*4+0x10]
vmovups XMMWORD PTR [rdi+rcx*4],xmm0
vmovups xmm0,XMMWORD PTR [rsi+rcx*4+0x20]
vaddps xmm0,xmm0,XMMWORD PTR [rdx+rcx*4+0x20]
vmovups XMMWORD PTR [rdi+rcx*4+0x10],xmm1
vmovups xmm1,XMMWORD PTR [rsi+rcx*4+0x30]
vaddps xmm1,xmm1,XMMWORD PTR [rdx+rcx*4+0x30]
vmovups XMMWORD PTR [rdi+rcx*4+0x20],xmm0
vmovups XMMWORD PTR [rdi+rcx*4+0x30],xmm1

但是这段代码仍然是第一个版本的一半性能......

【讨论】：

DRAM 时钟速度是内存带宽的一个潜在瓶颈，但一个核心的有限内存并行性是另一个问题：在现代台式机 CPU 上，单个核心甚至无法使双通道内存饱和，而且它是在多核 Xeon 上更糟糕：潜在的聚合带宽更高，但每核最大带宽更低。 Why is Skylake so much better than Broadwell-E for single-threaded memory throughput? 无论如何，很好的答案，但你应该提到你正在显示（我认为）clang 输出（不是 gcc）。我可以知道，因为 gcc 默认情况下根本不会展开。