uint64_t 与 int64_t 的 sqrt答案

【问题标题】：sqrt of uint64_t vs. int64_tuint64_t 与 int64_t 的 sqrt
【发布时间】：2018-05-20 20:30:54
【问题描述】：

我注意到计算uint64_t 平方根的整数部分比int64_t 复杂得多。请问，有人对此有解释吗？为什么处理一个额外的位似乎要困难得多？

以下内容：

int64_t sqrt_int(int64_t a) {
    return sqrt(a);
}

使用 clang 5.0 和 -mfpmath=sse -msse3 -Wall -O3 编译到

sqrt_int(long):                           # @sqrt_int(long)
        cvtsi2sd        xmm0, rdi
        sqrtsd  xmm0, xmm0
        cvttsd2si       rax, xmm0
        ret

但是如下：

uint64_t sqrt_int(uint64_t a) {
    return sqrt(a);
}

编译为：

.LCPI0_0:
        .long   1127219200              # 0x43300000
        .long   1160773632              # 0x45300000
        .long   0                       # 0x0
        .long   0                       # 0x0
.LCPI0_1:
        .quad   4841369599423283200     # double 4503599627370496
        .quad   4985484787499139072     # double 1.9342813113834067E+25
.LCPI0_2:
        .quad   4890909195324358656     # double 9.2233720368547758E+18
sqrt_int(unsigned long):                           # @sqrt_int(unsigned long)
        movq    xmm0, rdi
        punpckldq       xmm0, xmmword ptr [rip + .LCPI0_0] # xmm0 = xmm0[0],mem[0],xmm0[1],mem[1]
        subpd   xmm0, xmmword ptr [rip + .LCPI0_1]
        haddpd  xmm0, xmm0
        sqrtsd  xmm0, xmm0
        movsd   xmm1, qword ptr [rip + .LCPI0_2] # xmm1 = mem[0],zero
        movapd  xmm2, xmm0
        subsd   xmm2, xmm1
        cvttsd2si       rax, xmm2
        movabs  rcx, -9223372036854775808
        xor     rcx, rax
        cvttsd2si       rax, xmm0
        ucomisd xmm0, xmm1
        cmovae  rax, rcx
        ret

【问题讨论】：

正如 Zalman 指出的那样，它将 arg 转换为 double（因为这就是您的 C 源代码所做的），因此严格来说它不是 sqrt 的整数部分。并非所有大于 2^53 的整数都可以精确表示为 double，因此这可能涉及舍入。 en.wikipedia.org/wiki/Double-precision_floating-point_format

标签： assembly compilation x86 sse square-root

【解决方案1】：

首先，您需要清楚这段代码将 64 位整数（有符号或无符号）转换为双精度浮点数，取平方根，然后将结果转换回有符号或无符号整数。

您的问题的答案是因为英特尔在您正在编译的指令集中提供了 64 位有符号整数到双精度浮点转换（反之亦然），但对于无符号情况没有这样做。他们在 AVX-512 中添加了无符号转换指令，但在此之前它并不存在。因此，对于带符号的情况，转换为双精度和转换回各是一条指令。对于无符号情况，编译器必须根据可用指令综合转换操作。

您可以在此处获取有关 SSE2/AVX/AVX-512 等的哪些版本中可用的指令的信息： https://software.intel.com/sites/landingpage/IntrinsicsGuide/

您可以在此处查看有关用于合成转换的方法的讨论： Are there unsigned equivalents of the x87 FILD and SSE CVTSI2SD instructions?

【讨论】：

【解决方案2】：

除了 Zalman 出色的answer：sqrt 的结果总是小于 INT64_MAX，因为sqrt 的输入在uint64_t 范围内。因此，单个cvttsd2si 足以将双精度转换回uint64_t。换句话说：对于所有输入值a函数

uint64_t sqrt_int(uint64_t a) {
    return sqrt(a);
}

与修改后的代码具有完全相同的行为

uint64_t sqrt_int(uint64_t a) {
    return (int64_t)sqrt(a);
}

后一个函数编译成

.LCPI0_0:
  .long 1127219200 # 0x43300000
  .long 1160773632 # 0x45300000
  .long 0 # 0x0
  .long 0 # 0x0
.LCPI0_1:
  .quad 4841369599423283200 # double 4503599627370496
  .quad 4985484787499139072 # double 1.9342813113834067E+25
sqrt_int(unsigned long): # @sqrt_int(unsigned long)
  movq xmm0, rdi
  punpckldq xmm0, xmmword ptr [rip + .LCPI0_0] # xmm0 = xmm0[0],mem[0],xmm0[1],mem[1]
  subpd xmm0, xmmword ptr [rip + .LCPI0_1]
  haddpd xmm0, xmm0
  sqrtsd xmm0, xmm0
  cvttsd2si rax, xmm0
  ret

这比原始代码的指令少得多。

【讨论】：