提问人：Z boson 提问时间：6/1/2015 最后编辑：Peter CordesZ boson 更新时间：7/23/2021 访问量：1905
针对快速乘法但慢速加法进行优化：FMA 和 doubledouble
Optimize for fast multiplication but slow addition: FMA and doubledouble
问：
当我第一次获得 Haswell 处理器时，我尝试实现 FMA 来确定曼德布洛特集。主要算法是这样的：
intn = 0;
for(int32_t i=0; i<maxiter; i++) {
floatn x2 = square(x), y2 = square(y); //square(x) = x*x
floatn r2 = x2 + y2;
booln mask = r2<cut; //booln is in the float domain non integer domain
if(!horizontal_or(mask)) break; //_mm256_testz_pd(mask)
n -= mask
floatn t = x*y; mul2(t); //mul2(t): t*=2
x = x2 - y2 + cx;
y = t + cy;
}
这决定了像素是否在曼德布洛特集中。因此，对于双浮点，它运行在 4 个像素 （， ） 上。这需要 4 个 SIMD 浮点乘法和 4 个 SIMD 浮点加法。nfloatn = __m256dintn = __m256i
然后我修改了它以像这样与 FMA 一起使用
intn n = 0;
for(int32_t i=0; i<maxiter; i++) {
floatn r2 = mul_add(x,x,y*y);
booln mask = r2<cut;
if(!horizontal_or(mask)) break;
add_mask(n,mask);
floatn t = x*y;
x = mul_sub(x,x, mul_sub(y,y,cx));
y = mul_add(2.0f,t,cy);
}
mul_add呼叫和mul_sub呼叫。此方法使用 4 个 FMA SIMD 运算和两个 SIMD 乘法，与没有 FMA 相比，算术运算少了两个。此外，FMA 和乘法可以使用两个端口，并且只能添加一个端口。_mm256_fmad_pd_mm256_fmsub_pd
为了使我的测试不那么有偏差，我放大了一个完全在曼德布洛特集中的区域，所以所有的值都是 。在这种情况下，使用 FMA 的方法速度提高了约 27%。这当然是一个进步，但从 SSE 到 AVX 使我的性能翻了一番，所以我希望 FMA 能再增加两倍。maxiter
但后来我找到了关于 FMA 的答案，它说
融合乘加指令的重要方面是中间结果的（几乎）无限精度。这有助于提高性能，但不会那么大，因为两个操作被编码在一条指令中——它有助于提高性能，因为中间结果的几乎无限的精度有时很重要，并且当这种精度水平确实是程序员所追求的时，通过普通的乘法和加法来恢复非常昂贵。
后面给出了一个 double*double 到 double-double 乘法的例子
high = a * b; /* double-precision approximation of the real product */
low = fma(a, b, -high); /* remainder of the real product */
由此，我得出结论，我实现 FMA 的方式不是最佳的，因此我决定实施 SIMD double-double。我基于论文 Extended-Precision Floating-Point Numbers for GPU Computation 实现了 double-double。这张纸是用于双浮点的，所以我把它修改为双双。此外，我没有在 SIMD 寄存器中打包一个双双值，而是将 4 个双双值打包到一个 AVX 高电平寄存器和一个 AVX 低电平寄存器中。
对于曼德布洛特集，我真正需要的是双倍乘法和加法。在那篇论文中，这些是 和 函数。
下图显示了软件 FMA（左）和硬件 FMA（右）函数的程序集。这清楚地表明，硬件FMA是双倍乘法的一大改进。df64_adddf64_multdf64_mult
那么硬件FMA在双双曼德布洛特集计算中是如何表现的呢？答案是，这仅比软件 FMA 快 15% 左右。这比我所希望的要少得多。双倍曼德布洛特计算需要 4 次双倍加法和 4 次双倍乘法 （、、 和 ）。但是，2*（x*y） 乘法对于双倍是微不足道的，因此这种乘法在成本中可以忽略不计。因此，我认为使用硬件 FMA 的改进如此之小的原因是计算以缓慢的双双加法为主（见下面的组装）。x*xy*yx*y2*(x*y)
过去，乘法比加法慢（程序员使用了几种技巧来避免乘法），但对于 Haswell，情况似乎正好相反。不仅是由于 FMA，还因为乘法可以使用两个端口，但只能添加一个端口。
所以我的问题（最后）是：
- 与乘法相比，当加法较慢时，如何优化？
- 有没有一种代数方法可以改变我的算法以使用更多的乘法
和更少的添加？我知道有相反的方法，例如 多用两个加法少乘法。
(x+y)*(x+y) - (x*x+y*y) = 2*x*y - 有没有办法简化df64_add功能（例如使用FMA）？
如果有人想知道，double-double 方法比 double 慢十倍左右。这还不错，我认为如果有硬件四精度类型，它可能至少是双精度的两倍，所以我的软件方法比我对硬件的期望慢五倍（如果存在的话）。
df64_add集会
vmovapd 8(%rsp), %ymm0
movq %rdi, %rax
vmovapd 72(%rsp), %ymm1
vmovapd 40(%rsp), %ymm3
vaddpd %ymm1, %ymm0, %ymm4
vmovapd 104(%rsp), %ymm5
vsubpd %ymm0, %ymm4, %ymm2
vsubpd %ymm2, %ymm1, %ymm1
vsubpd %ymm2, %ymm4, %ymm2
vsubpd %ymm2, %ymm0, %ymm0
vaddpd %ymm1, %ymm0, %ymm2
vaddpd %ymm5, %ymm3, %ymm1
vsubpd %ymm3, %ymm1, %ymm6
vsubpd %ymm6, %ymm5, %ymm5
vsubpd %ymm6, %ymm1, %ymm6
vaddpd %ymm1, %ymm2, %ymm1
vsubpd %ymm6, %ymm3, %ymm3
vaddpd %ymm1, %ymm4, %ymm2
vaddpd %ymm5, %ymm3, %ymm3
vsubpd %ymm4, %ymm2, %ymm4
vsubpd %ymm4, %ymm1, %ymm1
vaddpd %ymm3, %ymm1, %ymm0
vaddpd %ymm0, %ymm2, %ymm1
vsubpd %ymm2, %ymm1, %ymm2
vmovapd %ymm1, (%rdi)
vsubpd %ymm2, %ymm0, %ymm0
vmovapd %ymm0, 32(%rdi)
vzeroupper
ret
答：
5赞
Z boson
6/4/2015
#1
为了回答我的第三个问题，我找到了一个更快的双倍加法解决方案。我在论文《实现》中找到了另一种定义 图形上的浮点-浮点运算符 硬件。
Theorem 5 (Add22 theorem) Let be ah+al and bh+bl the float-float arguments of the following
algorithm:
Add22 (ah ,al ,bh ,bl)
1 r = ah ⊕ bh
2 if | ah | ≥ | bh | then
3 s = ((( ah ⊖ r ) ⊕ bh ) ⊕ b l ) ⊕ a l
4 e l s e
5 s = ((( bh ⊖ r ) ⊕ ah ) ⊕ a l ) ⊕ b l
6 ( rh , r l ) = add12 ( r , s )
7 return (rh , r l)
这是我如何实现这个（伪代码）：
static inline doubledoublen add22(doubledoublen const &a, doubledouble const &b) {
doublen aa,ab,ah,bh,al,bl;
booln mask;
aa = abs(a.hi); //_mm256_and_pd
ab = abs(b.hi);
mask = aa >= ab; //_mm256_cmple_pd
// z = select(cut,x,y) is a SIMD version of z = cut ? x : y;
ah = select(mask,a.hi,b.hi); //_mm256_blendv_pd
bh = select(mask,b.hi,a.hi);
al = select(mask,a.lo,b.lo);
bl = select(mask,b.lo,a.lo);
doublen r, s;
r = ah + bh;
s = (((ah - r) + bh) + bl ) + al;
return two_sum(r,s);
}
Add22 的这个定义使用 11 个加法而不是 20 个加法，但它需要一些额外的代码来确定是否 .这里讨论如何实现 SIMD minmag 和 maxmag 函数。幸运的是，大多数附加代码不使用端口 1。现在只有 12 条指令进入端口 1，而不是 20 条。|ah| >= |bh|
这是新 Add22 的吞吐量分析表 IACA
Throughput Analysis Report
--------------------------
Block Throughput: 12.05 Cycles Throughput Bottleneck: Port1
Port Binding In Cycles Per Iteration:
---------------------------------------------------------------------------------------
| Port | 0 - DV | 1 | 2 - D | 3 - D | 4 | 5 | 6 | 7 |
---------------------------------------------------------------------------------------
| Cycles | 0.0 0.0 | 12.0 | 2.5 2.5 | 2.5 2.5 | 2.0 | 10.0 | 0.0 | 2.0 |
---------------------------------------------------------------------------------------
| Num Of | Ports pressure in cycles | |
| Uops | 0 - DV | 1 | 2 - D | 3 - D | 4 | 5 | 6 | 7 | |
---------------------------------------------------------------------------------
| 1 | | | 0.5 0.5 | 0.5 0.5 | | | | | | vmovapd ymm3, ymmword ptr [rip]
| 1 | | | 0.5 0.5 | 0.5 0.5 | | | | | | vmovapd ymm0, ymmword ptr [rdx]
| 1 | | | 0.5 0.5 | 0.5 0.5 | | | | | | vmovapd ymm4, ymmword ptr [rsi]
| 1 | | | | | | 1.0 | | | | vandpd ymm2, ymm4, ymm3
| 1 | | | | | | 1.0 | | | | vandpd ymm3, ymm0, ymm3
| 1 | | 1.0 | | | | | | | CP | vcmppd ymm2, ymm3, ymm2, 0x2
| 1 | | | 0.5 0.5 | 0.5 0.5 | | | | | | vmovapd ymm3, ymmword ptr [rsi+0x20]
| 2 | | | | | | 2.0 | | | | vblendvpd ymm1, ymm0, ymm4, ymm2
| 2 | | | | | | 2.0 | | | | vblendvpd ymm4, ymm4, ymm0, ymm2
| 1 | | | 0.5 0.5 | 0.5 0.5 | | | | | | vmovapd ymm0, ymmword ptr [rdx+0x20]
| 2 | | | | | | 2.0 | | | | vblendvpd ymm5, ymm0, ymm3, ymm2
| 2 | | | | | | 2.0 | | | | vblendvpd ymm0, ymm3, ymm0, ymm2
| 1 | | 1.0 | | | | | | | CP | vaddpd ymm3, ymm1, ymm4
| 1 | | 1.0 | | | | | | | CP | vsubpd ymm2, ymm1, ymm3
| 1 | | 1.0 | | | | | | | CP | vaddpd ymm1, ymm2, ymm4
| 1 | | 1.0 | | | | | | | CP | vaddpd ymm1, ymm1, ymm0
| 1 | | 1.0 | | | | | | | CP | vaddpd ymm0, ymm1, ymm5
| 1 | | 1.0 | | | | | | | CP | vaddpd ymm2, ymm3, ymm0
| 1 | | 1.0 | | | | | | | CP | vsubpd ymm1, ymm2, ymm3
| 2^ | | | | | 1.0 | | | 1.0 | | vmovapd ymmword ptr [rdi], ymm2
| 1 | | 1.0 | | | | | | | CP | vsubpd ymm0, ymm0, ymm1
| 1 | | 1.0 | | | | | | | CP | vsubpd ymm1, ymm2, ymm1
| 1 | | 1.0 | | | | | | | CP | vsubpd ymm3, ymm3, ymm1
| 1 | | 1.0 | | | | | | | CP | vaddpd ymm0, ymm3, ymm0
| 2^ | | | | | 1.0 | | | 1.0 | | vmovapd ymmword ptr [rdi+0x20], ymm0
这是旧的吞吐量分析
Throughput Analysis Report
--------------------------
Block Throughput: 20.00 Cycles Throughput Bottleneck: Port1
Port Binding In Cycles Per Iteration:
---------------------------------------------------------------------------------------
| Port | 0 - DV | 1 | 2 - D | 3 - D | 4 | 5 | 6 | 7 |
---------------------------------------------------------------------------------------
| Cycles | 0.0 0.0 | 20.0 | 2.0 2.0 | 2.0 2.0 | 2.0 | 0.0 | 0.0 | 2.0 |
---------------------------------------------------------------------------------------
| Num Of | Ports pressure in cycles | |
| Uops | 0 - DV | 1 | 2 - D | 3 - D | 4 | 5 | 6 | 7 | |
---------------------------------------------------------------------------------
| 1 | | | 1.0 1.0 | | | | | | | vmovapd ymm0, ymmword ptr [rsi]
| 1 | | | | 1.0 1.0 | | | | | | vmovapd ymm1, ymmword ptr [rdx]
| 1 | | | 1.0 1.0 | | | | | | | vmovapd ymm3, ymmword ptr [rsi+0x20]
| 1 | | 1.0 | | | | | | | CP | vaddpd ymm4, ymm0, ymm1
| 1 | | | | 1.0 1.0 | | | | | | vmovapd ymm5, ymmword ptr [rdx+0x20]
| 1 | | 1.0 | | | | | | | CP | vsubpd ymm2, ymm4, ymm0
| 1 | | 1.0 | | | | | | | CP | vsubpd ymm1, ymm1, ymm2
| 1 | | 1.0 | | | | | | | CP | vsubpd ymm2, ymm4, ymm2
| 1 | | 1.0 | | | | | | | CP | vsubpd ymm0, ymm0, ymm2
| 1 | | 1.0 | | | | | | | CP | vaddpd ymm2, ymm0, ymm1
| 1 | | 1.0 | | | | | | | CP | vaddpd ymm1, ymm3, ymm5
| 1 | | 1.0 | | | | | | | CP | vsubpd ymm6, ymm1, ymm3
| 1 | | 1.0 | | | | | | | CP | vsubpd ymm5, ymm5, ymm6
| 1 | | 1.0 | | | | | | | CP | vsubpd ymm6, ymm1, ymm6
| 1 | | 1.0 | | | | | | | CP | vaddpd ymm1, ymm2, ymm1
| 1 | | 1.0 | | | | | | | CP | vsubpd ymm3, ymm3, ymm6
| 1 | | 1.0 | | | | | | | CP | vaddpd ymm2, ymm4, ymm1
| 1 | | 1.0 | | | | | | | CP | vaddpd ymm3, ymm3, ymm5
| 1 | | 1.0 | | | | | | | CP | vsubpd ymm4, ymm2, ymm4
| 1 | | 1.0 | | | | | | | CP | vsubpd ymm1, ymm1, ymm4
| 1 | | 1.0 | | | | | | | CP | vaddpd ymm0, ymm1, ymm3
| 1 | | 1.0 | | | | | | | CP | vaddpd ymm1, ymm2, ymm0
| 1 | | 1.0 | | | | | | | CP | vsubpd ymm2, ymm1, ymm2
| 2^ | | | | | 1.0 | | | 1.0 | | vmovapd ymmword ptr [rdi], ymm1
| 1 | | 1.0 | | | | | | | CP | vsubpd ymm0, ymm0, ymm2
| 2^ | | | | | 1.0 | | | 1.0 | | vmovapd ymmword ptr [rdi+0x20], ymm0
更好的解决方案是，如果除了 FMA 之外还有三个操作数单舍入模式指令。在我看来，应该有单一的舍入模式说明
a + b + c
a * b + c //FMA - this is the only one in x86 so far
a * b * c
评论
1赞
tmyklebu
6/4/2015
什么能让你得到？三和使双倍加法变得微不足道，FMA使双倍乘法变得微不足道，但似乎没有明显的目的。a*b*ca*b*c
0赞
Z boson
6/4/2015
@tmyklebu，我不知道。也许你是对的。对于双倍，我同意你的看法，这似乎没有意义。我想知道这三个操作数一个舍入模式指令还有什么其他应用？我只用过双双。a*b*c
0赞
Z boson
6/12/2015
@tmyklebu，实际上想想看，唯一必要的是和.a + b - ca*b - c
1赞
Pierre
6/27/2018
#2
为了加快算法速度，我使用基于 2 fma、1 mul 和 2 add 的简化版本。我以这种方式处理 8 次迭代。然后计算逃逸半径，并在必要时回滚最后 8 次迭代。
编译器很好地展开了使用 x86 内部函数编写的以下临界循环 X = X^2 + C，在展开后，您会发现 2 个 FMA 操作彼此之间并没有严重依赖。
// IACA_START;
for (j = 0; j < 8; j++) {
Xrm = _mm256_mul_ps(Xre, Xim);
Xtt = _mm256_fmsub_ps(Xim, Xim, Cre);
Xrm = _mm256_add_ps(Xrm, Xrm);
Xim = _mm256_add_ps(Cim, Xrm);
Xre = _mm256_fmsub_ps(Xre, Xre, Xtt);
} // for
// IACA_END;
然后我计算逃逸半径（|X|< 阈值），这需要每 8 次迭代进行一次 FMA 和另一次乘法。
cmp = _mm256_mul_ps(Xre, Xre);
cmp = _mm256_fmadd_ps(Xim, Xim, cmp);
cmp = _mm256_cmp_ps(cmp, vec_threshold, _CMP_LE_OS);
if (_mm256_testc_si256((__m256i) cmp, vec_one)) {
i += 8;
continue;
}
你提到“加法很慢”，这并不完全正确，但你是对的，随着时间的推移，在最近的架构上，乘法吞吐量越来越高。
乘法延迟和依赖关系是关键。FMA 的吞吐量为 1 个周期，延迟为 5 个周期。独立 FMA 指令的执行可以重叠。
基于乘法结果的加法将获得完整的延迟命中。
因此，您必须通过执行“代码拼接”来打破这些直接依赖关系，并在同一循环中计算 2 个点，并且在与 IACA 检查将要发生的事情之前交错代码。以下代码有 2 组变量（后缀为 0 和 1 表示 X0=X0^2+C0、X1=X1^2+C1）并开始填充 FMA 孔
for (j = 0; j < 8; j++) {
Xrm0 = _mm256_mul_ps(Xre0, Xim0);
Xrm1 = _mm256_mul_ps(Xre1, Xim1);
Xtt0 = _mm256_fmsub_ps(Xim0, Xim0, Cre);
Xtt1 = _mm256_fmsub_ps(Xim1, Xim1, Cre);
Xrm0 = _mm256_add_ps(Xrm0, Xrm0);
Xrm1 = _mm256_add_ps(Xrm1, Xrm1);
Xim0 = _mm256_add_ps(Cim0, Xrm0);
Xim1 = _mm256_add_ps(Cim1, Xrm1);
Xre0 = _mm256_fmsub_ps(Xre0, Xre0, Xtt0);
Xre1 = _mm256_fmsub_ps(Xre1, Xre1, Xtt1);
} // for
总而言之，
- 您可以将关键循环中的指令数量减半
- 您可以添加更多独立指令，并利用乘法和融合乘法加法的高吞吐量与低延迟。
评论
0赞
Peter Cordes
6/27/2018
你所说的“代码拼接”看起来像软件流水线。但是，是的，保持更多数据在传输中以隐藏 FP 延迟会暴露更多的指令级并行性。展开更多点是将原始问题中的数据并行性转化为 SIMD + ILP 的好方法。
0赞
Peter Cordes
6/27/2018
您的 FMA 延迟/吞吐量数字似乎来自 Ryzen 或 Bulldozer 系列。但是，您的意思是“基于乘法结果的加法会获得完整的延迟效果”？任何依赖指令都会受到完全延迟的影响;FMA->FMA 具有与 MUL->ADD 相同的延迟（即使在 Bulldozer 上，FMA 单元内有特殊的快进功能，使其为 5c 而不是 6c。 agner.org/optimize）。无论如何，所有具有 FMA 的 Intel CPU 的 FMA 吞吐量为 0.5c。但是 ADD 的吞吐量为 1c，在 Skylake 之前的 CPU 上（它在 FMA 单元上执行所有操作，将其降低到 4c lat）。
0赞
Pierre
6/27/2018
0.5c 基于带有 2 个 FMA 引擎的单元。很难为 haswell、skylake、coffeelake、cubilake 的每个变体“填充管道”......我使用简单的阈值“理想的fma数”=“延迟”/“吞吐量”。我知道我建议的代码不是最佳的，因为超线程仍然可以填补漏洞。如果编译器可以避免将实时寄存器溢出到内存中（使用 skylake 的 avx512vl 变体......
0赞
Peter Cordes
6/28/2018
您使用的是 256 位 FMA。都支持 FMA 的 Intel CPU 有两个 256 位 FMA 单元。您正在考虑 Skylake-AVX512，其中一些只有一个 512 位 FMA 单元。（当 512 位 uops 在飞行中时，SKX 会关闭端口 1 上的矢量 ALU。端口 0/端口 1 上的两个 256 位 FMA 单元作为单个 512 位 FMA 单元工作。某些型号在端口 5 上有一个额外的 512 位 FMA 单元，该单元仅针对 512 位 uops 激活。
1赞
Z boson
6/28/2018
@Pierre 我问题的重点是关于使用以提高精度。我已经知道如何优化.放大过去后，您的方法将失去精度。但是我找到了一种更好的曼德布洛特集合方法，使用扰动理论将限制因子从尾数移动到指数。所以现在我可以去哪个比哪个好得多，而不会改变指数的精度。无论如何，据我所知，您的答案根本没有回答我的问题。double-doubledouble10^-710^-308double-double
1赞
Pierre
7/3/2018
#3
您提到以下代码：
vsubpd %ymm0, %ymm4, %ymm2
vsubpd %ymm2, %ymm1, %ymm1 <-- immediate dependency ymm2
vsubpd %ymm2, %ymm4, %ymm2
vsubpd %ymm2, %ymm0, %ymm0 <-- immediate dependency ymm2
vaddpd %ymm1, %ymm0, %ymm2 <-- immediate dependency ymm0
vaddpd %ymm5, %ymm3, %ymm1
vsubpd %ymm3, %ymm1, %ymm6 <-- immediate dependency ymm1
vsubpd %ymm6, %ymm5, %ymm5 <-- immediate dependency ymm6
vsubpd %ymm6, %ymm1, %ymm6 <-- dependency ymm1, ymm6
vaddpd %ymm1, %ymm2, %ymm1
vsubpd %ymm6, %ymm3, %ymm3 <-- dependency ymm6
vaddpd %ymm1, %ymm4, %ymm2
vaddpd %ymm5, %ymm3, %ymm3 <-- dependency ymm3
vsubpd %ymm4, %ymm2, %ymm4
vsubpd %ymm4, %ymm1, %ymm1 <-- immediate dependency ymm4
vaddpd %ymm3, %ymm1, %ymm0 <-- immediate dependency ymm1, ymm3
vaddpd %ymm0, %ymm2, %ymm1 <-- immediate dependency ymm0
vsubpd %ymm2, %ymm1, %ymm2 <-- immediate dependency ymm1
如果仔细检查，这些大多是依赖操作，并且不满足有关延迟/吞吐量效率的基本规则。大多数指令取决于前一个指令的结果，或之前的 2 个指令。此序列包含 30 个周期的关键路径（大约 9 或 10 条关于“3 个周期延迟”/“1 个周期吞吐量”的指令）。
您的 IACA 在关键路径中报告“CP”=> 条指令，评估成本为 20 个周期的吞吐量。您应该获得延迟报告，因为如果您对执行速度感兴趣，它很重要。
为了消除这个关键路径的成本，如果编译器无法做到这一点，你必须交错大约 20 个类似的指令（例如，因为你的双精度代码在一个单独的库中编译，没有 -flto 优化和 vzeroupper 在函数进入和退出时无处不在，矢量化器只适用于内联代码）。
一种可能性是并行运行 2 个计算（请参阅上一篇文章中关于代码拼接以改进流水线）
如果我假设你的双双代码看起来像这个“标准”实现
// (r,e) = x + y
#define two_sum(x, y, r, e)
do { double t; r = x + y; t = r - x; e = (x - (r - t)) + (y - t); } while (0)
#define two_difference(x, y, r, e) \
do { double t; r = x - y; t = r - x; e = (x - (r - t)) - (y + t); } while (0)
.....
然后，您必须考虑以下代码，其中指令以相当细的粒度交错。
// (r1, e1) = x1 + y1, (r2, e2) x2 + y2
#define two_sum(x1, y1, x2, y2, r1, e1, r2, e2)
do { double t1, t2 \
r1 = x1 + y1; r2 = x2 + y2; \
t1 = r1 - x1; t2 = r2 - x2; \
e1 = (x1 - (r1 - t1)) + (y1 - t1); e2 = (x2 - (r2 - t2)) + (y2 - t2); \
} while (0)
....
然后，这将创建如下代码（延迟报告中的关键路径大致相同，大约有 35 条指令）。有关运行时的详细信息，Out-Of-Order 执行应该会飞越它而不会停滞。
vsubsd %xmm2, %xmm0, %xmm8
vsubsd %xmm3, %xmm1, %xmm1
vaddsd %xmm4, %xmm4, %xmm4
vaddsd %xmm5, %xmm5, %xmm5
vsubsd %xmm0, %xmm8, %xmm9
vsubsd %xmm9, %xmm8, %xmm10
vaddsd %xmm2, %xmm9, %xmm2
vsubsd %xmm10, %xmm0, %xmm0
vsubsd %xmm2, %xmm0, %xmm11
vaddsd %xmm14, %xmm4, %xmm2
vaddsd %xmm11, %xmm1, %xmm12
vsubsd %xmm4, %xmm2, %xmm0
vaddsd %xmm12, %xmm8, %xmm13
vsubsd %xmm0, %xmm2, %xmm11
vsubsd %xmm0, %xmm14, %xmm1
vaddsd %xmm6, %xmm13, %xmm3
vsubsd %xmm8, %xmm13, %xmm8
vsubsd %xmm11, %xmm4, %xmm4
vsubsd %xmm13, %xmm3, %xmm15
vsubsd %xmm8, %xmm12, %xmm12
vaddsd %xmm1, %xmm4, %xmm14
vsubsd %xmm15, %xmm3, %xmm9
vsubsd %xmm15, %xmm6, %xmm6
vaddsd %xmm7, %xmm12, %xmm7
vsubsd %xmm9, %xmm13, %xmm10
vaddsd 16(%rsp), %xmm5, %xmm9
vaddsd %xmm6, %xmm10, %xmm15
vaddsd %xmm14, %xmm9, %xmm10
vaddsd %xmm15, %xmm7, %xmm13
vaddsd %xmm10, %xmm2, %xmm15
vaddsd %xmm13, %xmm3, %xmm6
vsubsd %xmm2, %xmm15, %xmm2
vsubsd %xmm3, %xmm6, %xmm3
vsubsd %xmm2, %xmm10, %xmm11
vsubsd %xmm3, %xmm13, %xmm0
总结：
内联双双源代码：由于 ABI 约束，编译器和矢量化器无法跨函数调用进行优化，并且由于担心混叠而无法跨内存访问进行优化。
拼接代码以平衡吞吐量和延迟，并最大限度地提高 CPU 端口使用率（并最大限度地提高每个周期的指令数），只要编译器不会将太多寄存器溢出到内存中。
您可以使用 perf 实用程序（软件包 linux-tools-generic 和 linux-cloud-tools-generic）跟踪优化影响，以获取执行的指令数和每个周期的指令数。
评论
0赞
Peter Cordes
7/3/2018
无序执行将隐藏依赖链，只要它们不是循环携带的，或者如果存在任何指令级并行性，则部分隐藏它们。如果循环携带的 dep 链的延迟是瓶颈，而不是前端或执行端口，IACA 报告会说“依赖”或“延迟”（我忘了是哪个）。现代 x86 上的指令调度通常不是很重要。
0赞
Pierre
7/4/2018
@PeterCordes 你是对的，但在曼德布洛特示例的上下文中，下一个循环的执行取决于当前循环的结果。我建议的解决方案是开始在如此简单的循环中交错 2 个独立的依赖链（并且每个依赖链都已经矢量化）。在 2 个超线程上运行相同的代码是填补延迟漏洞的另一种方法。
评论
x*x + y*ydoublex2hi + x2lo + y2hi + y2lo + cxx2hiy2hicx2*x*y+cyxyhi = x+y; lo = (x-hi)+ylo