高效计算三个无符号整数的平均值(无溢出)

nju*_*ffa 24 c algorithm bit-manipulation micro-optimization extended-precision

有一个现有的问题“3 个长整数的平均值”,它特别关注三个有符号整数的平均值的有效计算。

然而,无符号整数的使用允许额外的优化不适用于上一个问题中涵盖的场景。这个问题是关于三个无符号整数的平均值的有效计算,其中平均值向零舍入,即在数学术语中我想计算?(a + b + c) / 3 ?。

计算此平均值的一种直接方法是

 avg = a / 3 + b / 3 + c / 3 + (a % 3 + b % 3 + c % 3) / 3;
Run Code Online (Sandbox Code Playgroud)

首先,现代优化编译器会将除法转换为具有倒数加移位的乘法,并将模运算转换为反向乘法和减法,其中反向乘法可能使用许多体系结构上可用的scale_add习语,例如leax86_64的,addlsl #n在ARM,iscadd在NVIDIA GPU。

在尝试以适用于许多常见平台的通用方式优化上述内容时,我观察到整数运算的成本通常在逻辑关系中?(添加|)?转移规模_添加MUL。这里的成本是指所有延迟、吞吐量限制和功耗。当处理的整数类型比本地寄存器宽度更宽时,例如在uint64_t32 位处理器上处理数据时,任何此类差异都会变得更加明显。

因此,我的优化策略是尽量减少指令数量,并在可能的情况下用“廉价”操作替换“昂贵”操作,同时不增加寄存器压力并为宽无序处理器保留可利用的并行性。

第一个观察结果是,我们可以通过首先应用产生和值和进位值的 CSA(进位保存加法器)将三个操作数的和减少为两个操作数的和,其中进位值的权重是和的两倍价值。在大多数处理器上,基于软件的 CSA 的成本是五个逻辑s。一些处理器,如 NVIDIA GPU,有一条LOP3指令可以一举计算三个操作数的任意逻辑表达式,在这种情况下,CSA 压缩为两个LOP3s(注意:我还没有说服 CUDA 编译器发出这两个LOP3s;它目前生产四个LOP3s!)。

第二个观察是,因为我们正在计算除以 3 的模,所以我们不需要反向乘法来计算它。我们可以改用dividend % 3= ((dividend / 3) + dividend) & 3,减少模一再加一个合乎逻辑的,因为我们已经有分工的结果。这是通用算法的一个实例:股息 % (2 n -1) = ((股息 / (2 n -1) + 股息) & (2 n -1)。

最后,对于修正项中的除以 3,(a % 3 + b % 3 + c % 3) / 3我们不需要通用除以 3 的代码。由于被除数非常小,在 [0, 6] 中,我们可以简化x / 3(3 * x) / 8只需要一个scale_add加上一个shift

下面的代码显示了我目前正在进行的工作。使用编译器资源管理器检查为各种平台生成的代码显示了我期望的紧凑代码(使用 编译时-O3)。

但是,在使用 Intel 13.x 编译器对 Ivy Bridge x86_64 机器上的代码计时时,一个缺陷变得明显:uint64_t与简单版本相比,虽然我的代码改善了延迟(数据从 18 个周期缩短到 15 个周期),但吞吐量却变差了(从对于uint64_t数据,每 6.8 个周期产生一个结果到每 8.5 个周期产生一个结果)。更仔细地查看汇编代码,很明显为什么会这样:我基本上设法将代码从大致三路并行降低到大致两路并行。

是否有一种通用的优化技术,对通用处理器特别是 x86 和 ARM 以及 GPU 的所有版本都有好处,可以保留更多的并行性?或者,是否有一种优化技术可以进一步减少总操作数以弥补减少的并行度?校正项的计算(tail在下面的代码中)似乎是一个很好的目标。简化(carry_mod_3 + sum_mod_3) / 2看起来很诱人,但对于九种可能组合中的一种,却得出了错误的结果。

 avg = a / 3 + b / 3 + c / 3 + (a % 3 + b % 3 + c % 3) / 3;
Run Code Online (Sandbox Code Playgroud)

Dav*_*tat 15

让我把帽子扔进戒指。我想在这里不要做任何太棘手的事情。

#include <stdint.h>

uint64_t average_of_three(uint64_t a, uint64_t b, uint64_t c) {
  uint64_t hi = (a >> 32) + (b >> 32) + (c >> 32);
  uint64_t lo = hi + (a & 0xffffffff) + (b & 0xffffffff) + (c & 0xffffffff);
  return 0x55555555 * hi + lo / 3;
}
Run Code Online (Sandbox Code Playgroud)

在下面关于不同拆分的讨论之后,这里有一个以三个按位与为代价节省乘法的版本:

T hi = (a >> 2) + (b >> 2) + (c >> 2);
T lo = (a & 3) + (b & 3) + (c & 3);
avg = hi + (hi + lo) / 3;
Run Code Online (Sandbox Code Playgroud)

  • 我喜欢这种简单性;看来我的大脑陷入了太多的纠结。由于我(= OP)正在使用英特尔编译器,因此在延迟和吞吐量方面,该变体实际上比迄今为止尝试过的任何其他变体在 Ivy Bridge 上对我来说效果更好。除了 `x86_64` 上的 `gcc` 具有奇怪的截断习惯用法之外,还为我在 godbolt 上检查的所有编译器提供了更严格的代码。因为我需要 `uint32_t` 和 `uint64_t` 版本,所以我概括如下: `const int shift = ((sizeof(T) * CHAR_BIT) / 2); const T 掩码 = (((T)1) &lt;&lt; 移位) - ((T)1); const T 乘数 = (((T)(~(T)0)) / 3) &gt;&gt; 移位;` (3认同)
  • 对于 x86-64,Clang 在这方面比 GCC 做得更好。gcc 花费更多指令,并使用 `mov same,same` 截断为 32 位,击败 [mov-elimination](/sf/ask/3091853971/ -免费-为什么-我不能-完全复制-这个)。https://godbolt.org/z/ejYxvT (例如 `mov esi,esi` 将 ESI 零扩展为 RSI 仍然会花费后端 uop,而 `mov r8d, esi` 则不会,这对于 OP 来说可能很重要Ivy Bridge CPU(前端比后端宽)比 Haswell 及更高版本更佳。)如果在关键路径上,也可以解决所有地方的延迟问题。 (2认同)
  • @njuffa:对于 IvyBridge 及更高版本上的 32 位 x86,lo = 8、hi = 24 位的分割让我们可以使用“movzx dword, byte”,它允许对字节源 movzx 进行 mov 消除。(在 Intel 上,但不是 AMD。仅适用于字节源,不适用于字源。[与我的第一条评论相同的链接](/sf/ask/3091853971/ -why-cant-i-regenerate-this-at-all) 了解详细信息。)如果数学仍然有效,这将节省后端 uops,并可能减少关键路径延迟。 (2认同)

Fal*_*ner 6

我不确定它是否符合您的要求,但也许它可以只计算结果,然后从溢出中修复错误:

T average_of_3 (T a, T b, T c)
{
    T r = ((T) (a + b + c)) / 3;
    T o = (a > (T) ~b) + ((T) (a + b) > (T) (~c));
    if (o) r += ((T) 0x5555555555555555) << (o - 1);
    T rem = ((T) (a + b + c)) % 3;
    if (rem >= (3 - o)) ++r;
    return r;
}
Run Code Online (Sandbox Code Playgroud)

[编辑] 这是我能想到的最好的无分支比较的版本。在我的机器上,这个版本的吞吐量实际上比 njuffa 的代码略高。__builtin_add_overflow(x, y, r)由 gcc 和 clang 支持,1如果和x + y溢出类型,则返回*r0否则返回,因此计算o相当于第一个版本中的可移植代码,但至少 gcc 使用内置生成更好的代码。

T average_of_3 (T a, T b, T c)
{
    T r = ((T) (a + b + c)) / 3;
    T rem = ((T) (a + b + c)) % 3;
    T dummy;
    T o = __builtin_add_overflow(a, b, &dummy) + __builtin_add_overflow((T) (a + b), c, &dummy);
    r += -((o - 1) & 0xaaaaaaaaaaaaaaab) ^ 0x5555555555555555;
    r += (rem + o + 1) >> 2;
    return r;
}
Run Code Online (Sandbox Code Playgroud)


Dav*_*tat 5

新答案,新思路。这是基于数学恒等式

floor((a+b+c)/3) = floor(x + (a+b+c - 3x)/3)
Run Code Online (Sandbox Code Playgroud)

这何时适用于机器整数和无符号除法?
当差异不换行时,即0 ? a+b+c - 3x ? T_MAX.

这个定义x很快,可以完成工作。

T avg3(T a, T b, T c) {
  T x = (a >> 2) + (b >> 2) + (c >> 2);
  return x + (a + b + c - 3 * x) / 3;
}
Run Code Online (Sandbox Code Playgroud)

奇怪的是,除非我这样做,否则 ICC 会插入一个额外的否定:

T avg3(T a, T b, T c) {
  T x = (a >> 2) + (b >> 2) + (c >> 2);
  return x + (a + b + c - (x + x * 2)) / 3;
}
Run Code Online (Sandbox Code Playgroud)

请注意,T必须至少有 5 位宽。

如果T是两个平台字长,那么可以通过省略.的低字来节省一些双字操作x

延迟更差但吞吐量可能略高的替代版本?

T lo = a + b;
T hi = lo < b;
lo += c;
hi += lo < c;
T x = (hi << (sizeof(T) * CHAR_BIT - 2)) + (lo >> 2);
avg = x + (T)(lo - 3 * x) / 3;
Run Code Online (Sandbox Code Playgroud)


Kev*_*inZ 5

我已经回答了您链接的问题,所以我只回答与此不同的部分:性能。

如果你真的关心性能,那么答案是:

( a + b + c ) / 3
Run Code Online (Sandbox Code Playgroud)

由于您关心性能,因此您应该对正在处理的数据大小有一个直觉。您不应该担心只有 3 个值的加法溢出(乘法是另一回事),因为如果您的数据已经足够大,可以使用所选数据类型的高位,那么无论如何您都有溢出的危险,应该使用更大的整数类型。如果你在 uint64_t 上溢出,那么你真的应该问自己为什么你需要准确地计算到 18 quintillion,也许考虑使用 float 或 double。

现在,说了这么多,我会给你我的实际答复:没关系。这个问题在现实生活中不会出现,当它出现时,性能并不重要。

如果您在 SIMD 中执行一百万次,这可能是一个真正的性能问题,因为在那里,您真的很想使用宽度较小的整数,并且您可能需要最后一点空间,但这不是您的问题。

  • 如果您偏执,您可以使用 x86 `add rdi, rsi` / `jc safe_version` / `add rdi, rdx` / `jc safe_version` / ... 来实现此功能以检测环绕​​。在现代 Intel CPU 上,这些 add/jcc 指令将宏融合到添加和分支微指令中。仍然比不检查而添加的效率低,但如果您预计永远不会或很少需要安全版本,那么这是一个很好的乐观快速路径。但是,是的,大多数时候这是正确的答案,而且这个问题是一个 XY 问题。如果您还没有 uint64_t,请在必要时将输入零扩展为 64 位。 (2认同)