首页 文章

为什么GCC在实现整数除法时使用乘以奇数的乘法？

提问于
浏览 次
174

我一直在阅读 div 和 mul 汇编操作，我决定通过在C中编写一个简单的程序来看它们的运行情况：

文件division.c

#include <stdlib.h>
#include <stdio.h>

int main()
{
    size_t i = 9;
    size_t j = i / 5;
    printf("%zu\n",j);
    return 0;
}

然后生成汇编语言代码：

gcc -S division.c -O0 -masm=intel

但是看着生成的 division.s 文件，它并没有't contain any div operations! Instead, it does some kind of black magic with bit shifting and magic numbers. Here'是一个计算 i/5 的代码片段：

mov     rax, QWORD PTR [rbp-16]   ; Move i (=9) to RAX
movabs  rdx, -3689348814741910323 ; Move some magic number to RDX (?)
mul     rdx                       ; Multiply 9 by magic number
mov     rax, rdx                  ; Take only the upper 64 bits of the result
shr     rax, 2                    ; Shift these bits 2 places to the right (?)
mov     QWORD PTR [rbp-8], rax    ; Magically, RAX contains 9/5=1 now, 
                                  ; so we can assign it to j

这里发生了什么？为什么海湾合作委员会根本不使用div？它如何产生这个神奇的数字以及为什么一切都有效？

c gcc assembly x86-64 integer-division

4 回答

  • 94

    整数除法是您可以在现代处理器上执行的最慢的算术运算之一，延迟可达数十个周期和吞吐量不佳 . （对于x86，请参阅Agner Fog's instruction tables and microarch guide） .

    如果您提前知道除数，则可以通过将其替换为具有相同效果的一组其他运算（乘法，加法和移位）来避免除法 . 即使需要进行多次操作，它通常仍然比整数除法本身快得多 .

    以这种方式实现C / 运算符而不是涉及 div 的多指令序列只是GCC 's default way of doing division by constants. It doesn' t需要跨操作进行优化，即使调试也不会改变任何内容 . （但是对于小代码大小使用 -Os 确实让GCC使用 div . ）使用乘法逆而不是除法就像使用 lea 而不是 mul 和 add

    因此，如果在编译时未知除数，则只会在输出中看到 div 或 idiv .

    有关编译器如何生成这些序列的信息，以及允许您自己生成它们的代码（除非您使用脑死亡编译器，否则几乎肯定是不必要的），请参阅libdivide .

    回复于
  • 8

    除以5与乘以1/5相同，再乘以4/5并向右移2位相同 . 有关的值是十六进制的 CCCCCCCCCCCCD ，如果放在一个十六进制点之后是4/5的二进制表示（即四分之四的二进制是 0.110011001100 重复 - 见下面的原因） . 我想你可以从这里拿走它！您可能想要查看fixed point arithmetic（但请注意它在最后四舍五入为整数 .

    至于为什么，乘法比除法快，当除数是固定的时，这是一条更快的路线 .

    有关其工作原理的详细说明，请参阅Reciprocal Multiplication, a tutorial，并按定点说明 . 它显示了查找倒数的算法如何工作，以及如何处理有符号的除法和模数 .

    让我们考虑一下为什么 0.CCCCCCCC... （十六进制）或 0.110011001100... 二进制为4/5 . 将二进制表示除以4（右移2位），我们将得到 0.001100110011... ，通过简单的检查可以添加原始的 0.111111111111... ，显然等于1，十进制中的 0.9999999... 等于1 . 因此，我们知道 x + x/4 = 1 ，所以 5x/4 = 1 ， x=4/5 . 然后将其表示为十六进制的 CCCCCCCCCCCCD 用于舍入（因为超出最后一个的二进制数字将是 1 ） .

    回复于
  • 47

    通常，乘法比除法快得多 . 因此，如果我们可以通过乘以倒数来逃避，我们可以通过常数显着加快除法

    皱纹是我们不能准确地表示倒数（除非除法是2的幂，但在这种情况下我们通常只能将除法转换为位移） . 因此，为了确保正确的答案，我们必须小心，我们的倒数中的错误不会导致我们的最终结果出错 .

    -3689348814741910323是0xCCCCCCCCCCCCCCCD，它是刚好超过4/5的值，以0.64的固定点表示 .

    当我们将64位整数乘以0.64定点数时，我们得到64.64的结果 . 我们将值截断为64位整数（有效地将其舍入为零），然后执行进一步的移位，除以4并再次截断 . 通过查看位级别，很明显我们可以将两个截断视为单个截断 .

    这显然给了我们至少近似除以5的近似值，但它是否给我们一个正确的答案正确舍入为零？

    为了得到准确的答案，错误需要足够小，不要将答案推到舍入边界 .

    除以5的确切答案将始终具有0,1 / 5,2 / 5,3 / 5或4/5的小数部分 . 因此，乘法和移位结果中的正误差小于1/5永远不会把结果推到一个四舍五入的边界 .

    我们常量中的误差是（1/5）* 2-64 . i的值小于264，因此乘法后的误差小于1/5 . 除以4后，误差小于（1/5）* 2-2 .

    （1/5）* 2-2 <1/5所以答案总是等于做一个精确的除法并向零舍入 .

    不幸的是，这对所有除数都不起作用 .

    如果我们试图将4/7表示为0.64固定点数，并且从零开始舍入，则最终得到（6/7）* 2-64的误差 . 乘以不到264的i值后，我们最终得到的误差不到6/7，除以4后，我们最终得到的误差略低于1.5 / 7，大于1/7 .

    因此，为了正确地实现7，我们需要乘以0.65的固定点数 . 我们可以通过乘以固定点数的低64位来实现，然后加上原始数字（这可能会溢出到进位）然后通过进位进行旋转 .

    回复于
  • 134

    这里是一个算法文档的链接，它生成我在Visual Studio中看到的值和代码（在大多数情况下），并且我假设仍然在GCC中用于将变量整数除以常数整数 .

    http://gmplib.org/~tege/divcnst-pldi94.pdf

    在本文中，uword有N位，udword有2N位，n =分子，d =分母= divisor，l最初设置为ceil（log2（d）），shpre是pre-shift（在乘法之前使用）= e = d中的尾随零位数，shpost是移位后（乘法后使用），prec是精度= N-e = N-shpre . 目标是使用预移位，乘法和后移位来优化n / d的计算 .

    向下滚动到图6.2，它定义了如何生成udword乘数（最大大小为N 1位），但没有清楚地解释该过程 . 我将在下面解释 .

    图4.2和图6.2显示了对于大多数除数，乘法器如何减小到N位或更小的乘数 . 公式4.5解释了如何导出用于处理图4.1和4.2中的N 1位乘法器的公式 .

    回到图6.2 . 只有当除数> 2 ^（N-1）（当ℓ== N）时，分子才能大于udword，在这种情况下，n / d的优化替换是比较（如果n> = d，q = 1 ，否则q = 0），因此不会生成乘数 . mlow和mhigh的初始值将是N 1位，并且可以使用两个udword / uword除法来产生每个N 1位值（mlow或mhigh） . 以64位模式使用X86为例：

    ; upper 8 bytes of numerator = 2^(ℓ) = (upper part of 2^(N+ℓ))
; lower 8 bytes of numerator for mlow  = 0
; lower 8 bytes of numerator for mhigh = 2^(N+ℓ-prec) = 2^(ℓ+shpre) = 2^(ℓ+e)
numerator dq    2 dup(?)        ;16 byte numerator
divisor   dq    1 dup(?)        ; 8 byte divisor
; ...
        mov     rcx,divisor
        mov     rdx,0
        mov     rax,numerator+8    ;upper 8 bytes of numerator
        div     rcx                ;after div, rax == 1
        mov     rax,numerator      ;lower 8 bytes of numerator
        div     rcx
        mov     rdx,1              ;rdx:rax = N+1 bit value = 65 bit value

    您可以使用GCC进行测试 . 你已经看到了如何处理j = i / 5 . 看看如何处理j = i / 7（应该是N 1位乘法器的情况） .

    回复于

相关问题