首页 > 解决方案 > 如何为浮点值实现 totalOrder 谓词?

问题描述

IEEE 754 规范在 §5.10 中定义了一个总顺序,我想在汇编中实现它。

Wikipedia 的描述来看,这听起来很像可以实现无分支或几乎无分支,但我还没有想出一个像样的方法;而且我在主要编程语言中找不到任何现有的符合规范的实现

当比较两个浮点数时,它充当 ≤ 操作,除了 totalOrder(−0, +0) ∧ ¬ totalOrder(+0, −0) 和同一浮点数的不同表示按它们的顺序排列指数乘以符号位。然后通过排序 -qNaN < -sNaN < numbers < +sNaN < +qNaN 将排序扩展到 NaN,同一类中的两个 NaN 之间的排序基于这些数据的整数有效负载乘以符号位。

首先检查 NaN 然后跳转到浮点比较或处理 NaN 情况是否有意义,或者将浮点值移动到整数寄存器并在那里执行所有操作是否更有意义?

(至少从阅读描述中,感觉规范作者努力允许使用整数指令直接实现。)

在 x86-64 处理器上实现浮点总顺序的“最佳”方法是什么?

标签: assemblyfloating-pointx86-64ieee-754micro-optimization

解决方案


如果您将 FP 位模式与符号/大小整数(包括-0 < +0和 NaN 位模式1 )进行比较,这一切都可以正常工作。这就是为什么像binary64 ( double)这样的 IEEE 格式使用有偏指数并将字段按该顺序排列的原因之一。nextafter(另一个是通过位模式++或在位模式上的方便实现--。)

这可以通过 2 的补码整数比较有效地实现:

  • 如果标志都被清除:非负数 Just Work
  • 如果只有一个设置了符号位:任何负数都小于任何非负数,包括-0.0 < +0.020x80000000 < 0x00000000的补码x <= yJust Works。
  • 如果两者都设置了符号位 ( (x&y)>>63):2 的补码x<y是符号/幅度 FP x>y。在 x86 asm 中,您可以避免移位而只查看 SF,或者使用 SIMD 元素的高位。

    在不搞砸==情况的情况下处理这个问题是很棘手的:你不能只x&y对结果进行 XOR 签名<;当他们比较相等时,这会翻转它。<=当两个输入都为负时,它会给你,但<对于其他情况。我不确定这是否可用于排序。


使用SSE4.2 pcmpgtq,您可以在其普通 XMM 寄存器中对双 FP 值进行操作,或者对 32 位浮点数的SSE2(x86-64 保证) pcmpgtd进行操作。(请注意,与以下pcmpgtq相比相对较慢pcmpgtd:更少的端口和更高的延迟 。https ://agner.org/optimize/。例如在 Skylake 上,1 p5 uop 具有 3c 延迟,而 pcmpgtd 和 pcmpeqq 是 p0/p1 的 1 uop 1 个周期延迟。)

我们不能仅使用一个pcmpgtq+ 符号修正来处理按位相等的情况。
x1 bitwise_eq x0无论输入是正还是负,pcmpgtq 结果都为 0. sign(x0&x1)无论您希望 0 或 1 表示、 还是>根据>=总顺序来翻转它,都会产生不一致的行为。但不幸的是,FP 比较的行为意味着我们必须对 FP-equal 进行特殊处理,而不仅仅是 FP-unordered。<<=-0.0 == +0.0

您不需要汇编,只需uint64_t在 C 中键入双关语,例如让编译器可能使用movq rax, xmm0,或将内在函数用于向量 regs。

但是,如果您使用的是 asm,则可以考虑在 ZF=1 上进行 FP 比较和分支,这将设置为 unordered 或 equal,然后才进行整数。如果您希望 NaN 和完全相等(包括+-0.0 == -+0.0)很少见,这可能会很好。请注意 ZF,CF,PF = 1,1,1 对于文档中ucomisd无序。所有 x86 FP 都以相同的方式比较设置标志,直接或通过fcom// fnstsw axlahf

例如,独立版本可能如下所示。(在内联时简化,例如直接使用jb而不是setb调用者分支):

totalOrder:   ; 0/1 integer in EAX = (xmm0 <= xmm1 totalOrder)
    xor      eax, eax
    ucomisd  xmm0, xmm1           ; ZF=0 implies PF=0 (ordered) so just check ZF
    jz    .compare_as_integer     ; unordered or FP-equal
     ; else CF accurately reflects the < or > (total) order of xmm0 vs. xmm1
    setb     al                 ; or branch with jb
    ret

;; SSE4.2, using AVX 3-operand versions.  Use movaps as needed for non-AVX
### Untested
        ; Used for unordered or FP-equal, including -0.0 == +0.0
        ; but also including -1.0 == -1.0 for example
 .compare_as_integer:          ; should work in general for any sign/magnitude integer
    vpcmpgtq xmm2, xmm1, xmm0     ; reversed order of comparison: x1>x0 == x0<x1
    vpand    xmm3, xmm1, xmm0     ; we only care about the MSB of the 64-bit integer
    vpxor    xmm2, xmm3           ; flip if x0 & x1 are negative

    vpcmpeqq xmm1, xmm0
    vpor     xmm2, xmm1
       ; top bits of XMM2 hold the boolean result for each SIMD element
       ; suitable for use with blendvpd

    vmovmskpd eax, xmm2           ; low bit of EAX = valid, high bit might be garbage
    and      eax, 1          ; optional depending on use-case
     ; EAX=1 if x0 bitwise_eq x1 or sign/magnitude x1 > x0
    ret

使用 AVX512VL,vpternlogq可以替换所有 3 种 AND/XOR/OR 操作;它可以实现 3 个输入的任意布尔函数。 (y_gt_x) ^ (x&y) | y_eq_x.

没有 SSE4.2,或者只是作为一个标量无分支策略,我想出了这个。(例如,如果值实际上在内存中,那么您可以只mov加载而不是movq从 XMM regs)。

;; works on its own, or as the fallback after ucomisd/jz
compare_as_integer:
    movq     rcx, xmm0
    movq     rsi, xmm1

    xor      eax, eax
    cmp      rcx, rsi
   ; je  bitwise equal special case would simplify the rest
    setl     al                 ; 2's complement x < y
    sete     dl
    and      rcx, rsi           ; maybe something with TEST / CMOVS?
    shr      rcx, 63
    xor      al, cl           ; flip the SETL result if both inputs were negative
    or       al, dl           ; always true on bitwise equal
    ret

EAX 的异或归零应该可以安全地读取 EAX,即使在 P6 系列上,在使用setl8 位xoror. (为什么 GCC 不使用部分寄存器?)。在大多数其他 CPU 上,这里唯一的缺点是对 RDX 的旧值的错误依赖,我之前没有破坏它sete dl。如果我先对 EDX 进行异或归零,我们就可以xor进入orEAX。

一个分支策略可以像这样工作:

;; probably slower unless data is predictable, e.g. mostly non-negative
compare_as_integer_branchy:
    movq     rcx, xmm0
    movq     rsi, xmm1

    xor      eax, eax       ; mov eax,1 with je to a ret wouldn't avoid partial-register stalls for setl al
    cmp      rcx, rsi
    je      .flip_result        ; return 1
    setl     al                 ; 2's complement x < y

    test     rcx, rsi
    js     .flip_result         ; if (x&y both negative)
    ret

.flip_result:         ; not bitwise EQ, and both inputs negative
    xor      al, 1
    ret

如果您愿意,可以混合和匹配其中的一部分;AND/SHR/XOR 可以沿非相等路径使用,而不是test+js.


如果在对结果进行分支的情况下将其内联,则可以将 common(?)-case(有限且不等于)分支放在特殊情况处理之前。但是特殊情况包括有序<,因此 ZF=1 上的希望可预测分支(包括 PF=1 无序情况)可能仍然是一个好主意。

    ucomisd  xmm1, xmm0
    ja       x1_gt_x0                ; CF==0 && ZF==0
    ; maybe unordered, maybe -0 vs +0, maybe just x1 < x0

脚注 1:NaN 编码作为总订单的一部分

FP 值(及其符号/幅度编码)围绕零对称。符号位始终是符号位,即使对于 NaN 也是如此,因此可以以相同的方式处理。

  • 最小的幅度当然是 +-0.0:所有指数和尾数位为零。
  • 次正规线有一个零指数场(最小值),这意味着尾数的前导零。显式部分非零。幅度与尾数成线性关系。(零实际上只是次正规的特例。)
  • 归一化的数字跨越 exponent = 1 到 exponent < max,这意味着尾数中的前导 1。一个指数内的最大值(所有尾数位设置)刚好低于 ++exponent;尾数 = 0 值:即从尾数到指数的进位增加 1 到下一个可表示的浮点数/双精度数
  • +- 无穷大有指数 = 全一,尾数 = 0
  • +- NaN 有指数 = 全一,尾数 = 非零
    • 在 x86 上,sNaN 清除了尾数的最高位。休息是在任何地方至少有 1 个设置位的有效负载(否则它是一个 Inf)。
    • 在 x86 上,qNaN 具有尾数集的最高位。休息是有效载荷

https://cwiki.apache.org/confluence/display/stdcxx/FloatingPoint(链接自Are the bit patterns of NaNs really hardware-dependent?)在其他几个 ISA 上显示了一些 sNaN 和 qNaN 编码。有些与 x86 不同,但 POWER 和 Alpha 具有为 qNaN 设置的尾数的 MSB,因此它们具有比任何 sNaN 更大的整数幅度。

PA-RISC 选择了相反的方式,因此在该(过时的)ISA 上实现全序谓词需要为 FP-compare 无序情况做额外的工作;如果它们中的任何一个是任何类型的 NaN,那么在进行整数比较之前,在这两个值中翻转该位可能会起作用。

(我之所以提到这一点,是因为相同的算法可以用于可能不是专门用于 x86 的高级语言。但是您可能只想保留它并始终以相同的方式处理相同的二进制位模式,即使这意味着 qNaN < sNaN 在某些平台上。您甚至只能通过手动编写位模式来获得 sNaN。)


PS:我知道“significand”在技术上更正确,但是“mantissa”的音节更少,我更喜欢它,并且在这种情况下很好理解。


推荐阅读