0

我需要计算两个向量的点积：uint64_t a[N], b[N];( N<=60) 包含 64 位无符号整数。正是这个循环：

unsigned __int128 ans = 0;
for(int i=0;i<N;++i)
    ans += (unsigned __int128) a[i] * (unsigned __int128) b[i];

ans将溢出，因此结果必须保存在一个宽整数中，如 256 位。但是由于 N<=60 我们可以将结果保持在 160 (64*2 + 32) 位整数中。

缓慢的解决方案：

  1. 手动处理溢出：
unsigned __int128 ans = 0;
uint64_t overflow = 0;
for(int i=0;i<N;++i){
    auto temp = ans;
    ans += (unsigned __int128) a[i] * (unsigned __int128) b[i];
    if(ans<temp) overflow++;
}

这很慢，因为添加if会使程序减慢〜 2.2 倍。

  1. 使用类库boost::multiprecision::uint256_t或 GMP。

可能快速的解决方案：如果我们在 64 位机器上使用汇编编程，那么可以使用 3 个 64 位寄存器执行加法，方法是使用后add跟从低位到高位。adcadc

但我不想求助于 ASM，因为它很难维护，而且不便携。

我的目标是让它快速且可维护。

编辑：彼得在他的评论中指出，clang 支持我adc在 gcc 仍然使用分支时使用的想法（手动处理溢出）。

4

1 回答 1

2

您绝对不需要 256 位累加器。整数之和最多N只能产生进位，因此 64 位进位计数器足以处理 2^64 个元素长的向量的点积。即 128 位乘积之和的总宽度只需 192 = 64*3 或 160 = 128 + 32 位。即一个额外的寄存器。 Noverflow

是的，为此最佳的 x86-64 asm 是mov-load from on array，mul或者mulx使用来自另一个数组的内存源，然后add+ adcinto ans，并adc reg, 0累积进位。

（可能有一些循环展开，可能有 2 个或更多累加器（3 个寄存器组）。如果为 Intel CPU 展开，可能会避免mul内存源的索引寻址模式，因此它可以微融合负载。不幸的是 GCC / clang不要创建相对于另一个数组索引一个数组的循环，以最小化循环开销（1 个指针增量），同时避免其中一个数组的索引寻址模式，因此您不会从编译器获得最佳 asm。但是叮当很好。）

clang9.0 为您的代码生成-O3 -march=broadwell -fno-unroll-loops. 即使使用 128 位整数， Clang 也能识别a<b出无符号进位的常用习惯用法，从而导致/ / （不幸的是，clang 的循环展开取消优化为; ; ;而不是，破坏了循环展开的好处！这是一个错过的优化应该报告的错误。）a+=badd reg,regadc reg,regadc reg,0movsetcmovzxaddadc reg,0

GCC 实际上在if()您的
overflow += sum<prod;. 但这对 GCC 的其他主要缺失优化没有帮助：实际上与cmp/sbbfor 进行128 位比较，sum < prod而不是仅使用 last 的 CF 输出adc。GCC 知道如何针对uint64_t（使用进位标志进行有效的 128 位加法）进行优化，但显然不适用于从__uint128_t.

可能无法通过更改源来手动控制 gcc 以避免错过优化，这是 GCC 开发人员必须在 GCC 中修复的错误。（因此，您应该在他们的 bugzilla https://gcc.gnu.org/bugzilla/ 上将其报告为错过的优化错误；包括指向此 Q&A 的链接）。尝试自己使用块进行全手动uint64_t会更糟：中间块具有进位和进位。这很难在 C 中正确编写，而且 GCC 不会将其优化为adc.

即使使用overflow += __builtin_add_overflow(sum, prod, &sum);也不能完全帮助，给我们相同的cmp/sbb/adcGCC 代码生成！ https://gcc.gnu.org/onlinedocs/gcc/Integer-Overflow-Builtins.html说“编译器将尽可能使用硬件指令来实现这些内置函数，例如在添加后溢出条件跳转，有条件的进位跳转等”，但显然它只是不知道如何处理 128 位整数。

#include <stdint.h>

static const int N = 2048;
using u128 = unsigned __int128;
u128 dp(uint64_t a[], uint64_t b[], uint64_t *overflow_out)
{
    u128 sum = 0;
    uint64_t overflow = 0;
    for(int i=0;i<N;++i){
        u128 prod = a[i] * (unsigned __int128) b[i];
        //sum += prod;
        //overflow += sum<prod;       // gcc uses adc after a 3x mov / cmp / sbb; clang is still good
        overflow += __builtin_add_overflow(sum, prod, &sum);  // gcc less bad, clang still good
    }

    *overflow_out = overflow;
    return sum;
}

Clang 很好地编译了这个（Godbolt）：

# clang9.0 -O3 -march=broadwell -fno-unroll-loops  -Wall -Wextra
     ... zero some regs, copy RDX to R8 (because MUL uses that implicitly)

.LBB0_1:                                # =>This Inner Loop Header: Depth=1
        mov     rax, qword ptr [rsi + 8*rcx]
        mul     qword ptr [rdi + 8*rcx]
        add     r10, rax
        adc     r9, rdx                 # sum += prod
        adc     r11, 0                  # overflow += carry_out
        inc     rcx
        cmp     rcx, 2048
        jne     .LBB0_1               # }while(i < N);

        mov     qword ptr [r8], r11   # store the carry count
        mov     rax, r10
        mov     rdx, r9               # return sum in RDX:RAX
        ret

请注意，您不需要ADOX / ADCX 或 MULX 来提高效率。乱序 exec 可以交错多个短 FLAGS 依赖链。

您可以将另外 3 个寄存器用于另一个 192 位累加器（求和和进位），但这可能是不必要的。

这看起来像前端的 9 微指令（假设mul不能保持索引寻址模式微融合（非层压），但这cmp/jcc会宏融合），因此它最多每 2.25 个时钟周期运行 1 次乘法。Haswell 和更早的版本adc reg,reg以 2 uop 运行，但 Broadwell 将 3 输入指令（FLAGS + 2 regs）改进为 1 uop。 adc reg,0在 SnB 系列上特例为 1 uop。

循环携带的依赖链每个只有 1 个循环：add进入 r10、adc进入 r9 和adc进入 r11。这些 ADC 指令的 FLAGS 输入是一个短的非循环携带依赖链的一部分，乱序执行将在早餐时吃掉它。

具有 5 宽管道的 Ice Lake 将运行得更快一些，例如每次迭代可能 1.8 个周期，假设它仍然 unlaminatemul的内存操作数。

Zen 有一个 5 条指令宽的管道，如果它包含任何多 uop 指令，则为 6 uop 宽。因此，它可能会在每次迭代 2 个周期时运行它，在其 2c 吞吐量方面遇到瓶颈以进行完全乘法。（正常imul r64, r64的非扩展乘法是 1/clock，Zen 上的延迟为 3c，与 Intel 相同。）

Zen 2 加速mul m64到mulx1/clock ( https://www.uops.info/html-tp/ZEN2/MUL_M64-Measurements.html )，并且可能是这个循环中最快的 CPU .

通过一些展开和索引一个数组相对于另一个数组的手动优化版本，每次乘法的成本可能接近 6 uop = mov-load (1 uop) + mul mem (2 uop) + add + 2x adc (3 uops），在 Broadwell 上大约 1.5 个周期/倍。

它仍然会成为前端的瓶颈，并且 2 uop（指针增量和 cmp/jcc）的最小循环开销意味着展开 4 可以给我们每 6.5 个周期而不是每 6 个周期的 4 次乘法，或 92%完全展开的理论最大值。（假设内存或次优调度没有停顿，或者至少无序 exec 足以吸收这一点并且不会停顿前端。后端应该能够跟上前端在这里，在 Haswell 和以后，所以 ROB 应该有空间来吸收一些小气泡，但可能不是 L3 或 TLB 未命中。）

我认为在这里使用 SIMD 没有任何好处。加法的最宽元素大小是 64 位，甚至 AVX512 也只有 64x64 => 128 位乘法。除非您可以转换为double在这种情况下 AVX512 可以运行得更快。

于 2020-01-03T13:20:25.803 回答