c - 我需要一个用于定点数学库的快速 96 位对 64 位特定除法算法

Question

我目前正在编写一个快速的 32.32 定点数学库。我成功地使加法、减法和乘法正确工作，但我在除法上遇到了困难。

给记不住的小伙伴们提个醒:32.32定点数是一个32位整数部分和32位小数部分的数。

我想出的最佳算法需要 96 位整数除法，这是编译器通常没有内置的东西。

无论如何，这里是:

G = 2^32

notation: x is the 64-bit fixed-point number, x1 is its low nibble and x2 is its high

G*(a/b) = ((a1 + a2*G) / (b1 + b2*G))*G      // Decompose this

G*(a/b) = (a1*G) / (b1*G + b2) + (a2*G*G) / (b1*G + b2)

如您所见，(a2*G*G) 保证大于常规的 64 位整数。如果我的编译器实际上支持 uint128_t，我将简单地执行以下操作:

((uint128_t)x << 32) / y)

嗯，他们不是，我需要一个解决方案。感谢您的帮助。

Answer 1

您可以将较大的除法分解为多个 block ，这些 block 用较少的位进行除法。正如另一位发帖人已经提到的，该算法可以在 Knuth 的 TAOCP 中找到。

不过，没必要买这本书!

hacker delight 网站上有一个用 C 语言实现该算法的代码。它被编写为仅使用 32 位算术进行 64 位无符号除法，因此您不能直接剪切'n'粘贴代码。要从 64 位到 128 位，您必须将所有类型、掩码和常量加宽两倍，例如short 变成 int，0xffff 变成 0xffffffffll 等。

经过这个简单的更改后，您应该能够进行 128 位除法。

代码镜像到GitHub , 但最初发布于 Hackersdelight.org (原始链接不再可用)。

由于您的最大值只需要 96 位，64 位除法之一将始终返回零，因此您甚至可以稍微简化代码。

哦 - 在我忘记这一点之前:该代码仅适用于无符号值。要将有符号除法转换为无符号除法，您可以执行以下操作(伪代码样式):

fixpoint Divide (fixpoint a, fixpoint b)
{
    // check if the integers are of different sign:
    fixpoint sign_difference = a ^ b; 
    
    // do unsigned division:
    fixpoint x = unsigned_divide (abs(a), abs(b));
    
    // if the signs have been different: negate the result.
    if (sign_difference < 0)
    {
       x = -x;
    }
    
    return x;
}

该网站本身也值得一试:http://www.hackersdelight.org/

顺便说一句 - 你正在做的很好的任务..你介意告诉我们你需要定点库吗？

---

顺便说一句 - 用于除法的普通移位和减法算法也可以工作。

如果您的目标是 x86，您可以使用 MMX 或 SSE 内在函数来实现它。该算法仅依赖原始操作，因此执行速度也相当快。