Lou*_*s15 5 performance sqrt visual-c++
为了学术目的，我在 C++ 中使用了基本的数学函数实现。今天，我对以下平方根代码进行了基准测试：
inline float sqrt_new(float n)
{
__asm {
fld n
fsqrt
}
}
我很惊讶地看到它始终比标准sqrt函数快（它需要标准函数执行时间的 85% 左右）。
我不太明白为什么，很想更好地理解它。下面我展示了我用来分析的完整代码（在 Visual Studio 2015 中，在发布模式下编译并打开所有优化）：
#include <iostream>
#include <random>
#include <chrono>
#define M 1000000
float ranfloats[M];
using namespace std;
inline float sqrt_new(float n)
{
__asm {
fld n
fsqrt
}
}
int main()
{
default_random_engine randomGenerator(time(0));
uniform_real_distribution<float> diceroll(0.0f , 1.0f);
chrono::high_resolution_clock::time_point start1, start2;
chrono::high_resolution_clock::time_point end1, end2;
float sqrt1 = 0;
float sqrt2 = 0;
for (int i = 0; i<M; i++) ranfloats[i] = diceroll(randomGenerator);
start1 = std::chrono::high_resolution_clock::now();
for (int i = 0; i<M; i++) sqrt1 += sqrt(ranfloats[i]);
end1 = std::chrono::high_resolution_clock::now();
start2 = std::chrono::high_resolution_clock::now();
for (int i = 0; i<M; i++) sqrt2 += sqrt_new(ranfloats[i]);
end2 = std::chrono::high_resolution_clock::now();
auto time1 = std::chrono::duration_cast<std::chrono::milliseconds>(end1 - start1).count();
auto time2 = std::chrono::duration_cast<std::chrono::milliseconds>(end2 - start2).count();
cout << "Time elapsed for SQRT1: " << time1 << " seconds" << endl;
cout << "Time elapsed for SQRT2: " << time2 << " seconds" << endl;
cout << "Average of Time for SQRT2 / Time for SQRT1: " << time2 / time1 << endl;
cout << "Equal to standard sqrt? " << (sqrt1 == sqrt2) << endl;
system("pause");
return 0;
}
编辑：我正在编辑问题以包括计算平方根的两个循环的反汇编代码，因为它们在 Visual Studio 2015 中出现。
一、拆解为for (int i = 0; i<M; i++) sqrt1 += sqrt(ranfloats[i]);：
00091194 0F 5A C0 cvtps2pd xmm0,xmm0
00091197 E8 F2 18 00 00 call __libm_sse2_sqrt_precise (092A8Eh)
0009119C F2 0F 5A C0 cvtsd2ss xmm0,xmm0
000911A0 83 C6 04 add esi,4
000911A3 F3 0F 58 44 24 4C addss xmm0,dword ptr [esp+4Ch]
000911A9 F3 0F 11 44 24 4C movss dword ptr [esp+4Ch],xmm0
000911AF 81 FE 90 5C 46 00 cmp esi,offset __dyn_tls_dtor_callback (0465C90h)
000911B5 7C D9 jl main+190h (091190h)
接下来，反汇编为for (int i = 0; i<M; i++) sqrt2 += sqrt_new(ranfloats[i]);：
00091290 F3 0F 10 00 movss xmm0,dword ptr [eax]
00091294 F3 0F 11 44 24 6C movss dword ptr [esp+6Ch],xmm0
0009129A D9 44 24 6C fld dword ptr [esp+6Ch]
0009129E D9 FA fsqrt
000912A0 D9 5C 24 6C fstp dword ptr [esp+6Ch]
000912A4 F3 0F 10 44 24 6C movss xmm0,dword ptr [esp+6Ch]
000912AA 83 C0 04 add eax,4
000912AD F3 0F 58 44 24 54 addss xmm0,dword ptr [esp+54h]
000912B3 F3 0F 11 44 24 54 movss dword ptr [esp+54h],xmm0
000912B9 ?? ?? ??
000912BA ?? ?? ??
000912BB ?? ?? ??
000912BC ?? ?? ??
000912BD ?? ?? ??
000912BE ?? ?? ??
000912BF ?? ?? ??
000912C0 ?? ?? ??
000912C1 ?? ?? ??
000912C2 ?? ?? ??
000912C3 ?? ?? ??
000912C4 ?? ?? ??
000912C5 ?? ?? ??
000912C6 ?? ?? ??
000912C7 ?? ?? ??
000912C8 ?? ?? ??
000912C9 ?? ?? ??
000912CA ?? ?? ??
000912CB ?? ?? ??
000912CC ?? ?? ??
000912CD ?? ?? ??
000912CE ?? ?? ??
000912CF ?? ?? ??
000912D0 ?? ?? ??
000912D1 ?? ?? ??
000912D2 ?? ?? ??
000912D3 ?? ?? ??
000912D4 ?? ?? ??
000912D5 ?? ?? ??
000912D6 ?? ?? ??
000912D7 ?? ?? ??
000912D8 ?? ?? ??
000912D9 ?? ?? ??
000912DA ?? ?? ??
000912DB ?? ?? ??
000912DC ?? ?? ??
000912DD ?? ?? ??
000912DE ?? ?? ??
您的两个循环都非常糟糕，除了 sqrt 函数调用或 FSQRT 指令之外，还有许多瓶颈。并且至少比最佳标量 SQRTSS（单精度）代码慢 2 倍。这可能比体面的 SSE2 矢量化循环可以实现的速度慢 8 倍。即使不重新排序任何数学运算，您也可以击败 SQRTSS 吞吐量。
https://gcc.gnu.org/wiki/DontUseInlineAsm 中的许多原因都适用于您的示例。编译器将无法通过您的函数传播常量，并且它不会知道结果总是非负的（如果它不是 NaN）。fabs()如果您稍后对数字进行平方，它也将无法将其优化为。
同样非常重要的是，您可以使用 SSE2 SQRTPS ( _mm_sqrt_ps())击败自动矢量化。使用内在函数的“无错误检查”标量 sqrt() 函数也存在该问题。IDK 是否有任何方法可以在没有 的情况下获得最佳结果/fp:fast，但我对此表示怀疑。（除了在汇编中编写整个循环，或者自己使用内在函数向量化整个循环）。
令人印象深刻的是，您的 Haswell CPU 设法尽可能快地运行函数调用循环，尽管 inline-asm 循环甚至可能不会使 FSQRT 吞吐量饱和。
出于某种原因，您的库函数调用正在调用double sqrt(double)，而不是 C++ 重载float sqrt(float)。这会导致转换为 double 并返回到 float。可能您需要#include <cmath>获取重载，或者您可以调用sqrtf(). Linux 上的 gcc 和 clang 使用您当前的代码调用 sqrtf()（不转换为 double 和 back），但也许它们的<random>标头恰好包含<cmath>，而 MSVC 没有。或者也许还有其他事情发生。
库函数调用循环将总和保存在内存中（而不是寄存器）。显然 32 位版本使用的调用约定__libm_sse2_sqrt_precise不保留任何 XMM 寄存器。Windows x64 ABI 确实保留了 XMM6-XMM15，但维基百科说这是新的，而 32 位 ABI 没有这样做。我假设如果有任何保留调用的 XMM 寄存器，MSVC 的优化器会利用它们。
无论如何，除了在每个独立的标量浮点数上调用 sqrt 的吞吐量瓶颈之外，对 sqrt1 的循环携带依赖是一个延迟瓶颈，包括存储转发往返：
000911A3 F3 0F 58 44 24 4C addss xmm0,dword ptr [esp+4Ch]
000911A9 F3 0F 11 44 24 4C movss dword ptr [esp+4Ch],xmm0
乱序执行让每次迭代的其余代码重叠，所以你只是吞吐量的瓶颈，但无论库 sqrt 函数多么有效，这个延迟瓶颈都会将循环限制为每 6 + 3 = 9 个循环一次迭代。（Haswell ADDSS 延迟 = 3，XMM 加载/存储的存储转发延迟 = 6 个周期。比整数寄存器的存储转发多 1 个周期。请参阅Agner Fog 的指令表。）
SQRTSD 每 8-14 个周期的吞吐量为 1，因此循环携带的依赖性不是 Haswell 上的限制瓶颈。
inline-asm 版本具有用于 sqrt 结果的存储/重新加载往返，但它不是循环携带的依赖链的一部分。 MSVC inline-asm 语法使得很难避免存储转发往返以将数据输入/输出内联 asm。但更糟糕的是，您在 x87 堆栈上生成结果，而编译器想要在 XMM 寄存器中进行 SSE 数学运算。
然后 MSVC 无缘无故地在脚下射击，将总和保存在内存中而不是 XMM 寄存器中。它查看 inline-asm 语句以查看它们影响哪些寄存器，因此 IDK 为什么它没有看到您的 inline-asm 语句不会破坏任何 XMM regs。
所以 MSVC在这里做的工作比必要的要糟糕得多：
00091290 movss xmm0,dword ptr [eax] # load from the array
00091294 movss dword ptr [esp+6Ch],xmm0 # store to the stack
0009129A fld dword ptr [esp+6Ch] # x87 load from stack
0009129E fsqrt
000912A0 fstp dword ptr [esp+6Ch] # x87 store to the stack
000912A4 movss xmm0,dword ptr [esp+6Ch] # SSE load from the stack (of sqrt(array[i]))
000912AA add eax,4
000912AD addss xmm0,dword ptr [esp+54h] # SSE load+add of the sum
000912B3 movss dword ptr [esp+54h],xmm0 # SSE store of the sum
因此它具有与函数调用循环相同的循环携带依赖链（ADDSS + 存储转发）。Haswell FSQRT 每 8-17 个周期有一个吞吐量，所以它可能仍然是瓶颈。（所有涉及数组值的存储/重新加载对于每次迭代都是独立的，无序执行可以重叠许多迭代以隐藏该延迟链。但是，它们会阻塞加载/存储执行单元，有时会延迟关键-path 加载/存储一个额外的周期。这称为资源冲突。）
如果没有/fp:fast，sqrtf()库函数必须设置errno结果是否为 NaN。这就是它不能内联到 SQRTSS 的原因。
如果您确实想自己实现一个无检查标量 sqrt 函数，您可以使用 Intel 内在语法来实现：
// DON'T USE THIS, it defeats auto-vectorization
static inline
float sqrt_scalar(float x) {
__m128 xvec = _mm_set_ss(x);
xvec = _mm_cvtss_f32(_mm_sqrt_ss(xvec));
}
这将编译为带有 gcc 和 clang（没有 -ffast-math）的近乎最佳的标量循环。在 Godbolt 编译器浏览器上查看：
# gcc6.2 -O3 for the sqrt_new loop using _mm_sqrt_ss. good scalar code, but don't optimize further.
.L10:
movss xmm0, DWORD PTR [r12]
add r12, 4
sqrtss xmm0, xmm0
addss xmm1, xmm0
cmp r12, rbx
jne .L10
这个循环应该只在 SQRTSS 吞吐量上出现瓶颈（在 Haswell 上每 7 个时钟一个，明显比 SQRTSD 或 FSQRT 快），并且没有资源冲突。但是，与即使不重新排序 FP 添加（因为FP add/mul 不是真正的关联）相比，它仍然是垃圾：智能编译器（或使用内在函数的程序员）将使用 SQRTPS 获得 4 个结果与 SQRTSS 产生的 1 结果相同的吞吐量。将 SQRT 结果的向量解包为 4 个标量，然后您可以保持完全相同的操作顺序和相同的中间结果舍入。我很失望 clang 和 gcc 没有这样做。
然而，gcc 和 clang 确实设法避免调用库函数。clang3.9（仅-O3使用 ）使用 SQRTSS 甚至不检查 NaN。我认为这是合法的，而不是编译器错误。也许它看到代码没有使用errno？
另一方面，gcc6.2 推测性地内联 sqrtf()，带有 SQRTSS 并检查输入以查看它是否需要调用库函数。
# gcc6.2's sqrt() loop, without -ffast-math.
# speculative inlining of SQRTSS with a check + fallback
# spills/reloads a lot of stuff to memory even when it skips the call :(
# xmm1 = 0.0 (gcc -fverbose-asm says it's holding sqrt2, which is zero-initialized, so I guess gcc decides to reuse that zero)
.L9:
movss xmm0, DWORD PTR [rbx]
sqrtss xmm5, xmm0
ucomiss xmm1, xmm0 # compare input against 0.0
movss DWORD PTR [rsp+8], xmm5
jbe .L8 # if(0.0 <= SQRTSS input || unordered(0.0, input)) { skip the function call; }
movss DWORD PTR [rsp+12], xmm1 # silly gcc, this store isn't needed. ucomiss doesn't modify xmm1
call sqrtf # called for negative inputs, but not for NaN.
movss xmm1, DWORD PTR [rsp+12]
.L8:
movss xmm4, DWORD PTR [rsp+4] # silly gcc always stores/reloads both, instead of putting the stores/reloads inside the block that the jbe skips
addss xmm4, DWORD PTR [rsp+8]
add rbx, 4
movss DWORD PTR [rsp+4], xmm4
cmp rbp, rbx
jne .L9
不幸的是，gcc 在这里用脚射击，就像 MSVC 对 inline-asm 所做的一样：有一个存储转发往返作为循环携带的依赖项。所有溢出/重新加载都可能在 JBE 跳过的块内。也许 gcc 的负面输入会很常见。
更糟糕的是，如果您确实使用了/fp:fastor -ffast-math，即使是像 clang 这样聪明的编译器也无法将您_mm_sqrt_ss的代码重写为 SQRTPS。Clang 通常不仅擅长将内在函数 1:1 映射到指令，而且如果您错过了组合事物的机会，它还会提出更优化的洗牌和混合。
所以在启用快速 FP 数学的情况下，使用_mm_sqrt_ss是一个很大的损失。clang 将sqrt()库函数调用版本编译为 RSQRTPS + 牛顿-拉夫森迭代。
另请注意，您的微基准测试代码对实现的延迟不敏感sqrt_new()，仅对吞吐量敏感。延迟在实际 FP 代码中通常很重要，而不仅仅是吞吐量。但在其他情况下，比如对许多数组元素独立做同样的事情，延迟并不重要，因为乱序执行可以通过在许多循环迭代中运行指令来很好地隐藏它。
正如我之前提到的，数据在进出 MSVC 风格的 inline-asm过程中产生的额外存储/重新加载往返的延迟在这里是一个严重的问题。当 MSVC 内联函数时，fld n不会直接来自数组。
顺便说一句，Skylake 的 SQRTPS/SS 吞吐量为每 3 个周期一个，但仍有 12 个周期的延迟。SQRTPD/SD 吞吐量 = 每 4-6 个周期一个，延迟 = 15-16 个周期。所以 FP 平方根在 Skylake 上比在 Haswell 上更加流水线化。这放大了基准 FP sqrt 延迟与吞吐量之间的差异。