是否可以在没有 bignums 的 JavaScript 中执行快速的 48 位无符号乘法？

Posted 2023-02-19

【中文标题】是否可以在没有 bignums 的 JavaScript 中执行快速的 48 位无符号乘法？

【英文标题】：Is it possible to perform fast 48-bit unsigned multiplication in JavaScript without bignums?

【发布时间】：2022-01-08 07:40:30

【问题描述】：

在 javascript 中，我们可以进行 48 位的加减除模：

In JavaScript, we can perform 48-bit addition, subtraction, division and modulus, using the native Number type:

function u48_add(a, b) 
  return (a + b) % Math.pow(2, 48);


function u48_sub(a, b) 
  return (a - b + Math.pow(2,48)) % Math.pow(2, 48);


function u48_div(a, b) 
  return Math.floor(a / b);


function u48_mod(a, b) 
  return a % b;

所有这些操作都有效，因为中间值不能通过Number.MAX_SAFE_INTEGER。但是，对于乘法，他们可以：

function u48_mul(a, b) 
  return (a * b) % Math.pow(2, 48);

所以u48_mul 可能会返回不正确的结果。一个解决方案是使用 BigInt：

function u48_mul(a, b) 
  return Number((BigInt(a) * BigInt(b)) % (2n ** 48n));

但是，在大多数浏览器中，它的速度要慢得多。有没有什么巧妙的技巧可以让我们在 JavaScript 中更快地执行 48 位无符号乘法？

【问题讨论】：

但是大于 2^53 的中间值是合理预期的吗？同样在您的 bigInt 示例中，我认为您会通过再次强制 Number 来降低精度

@johnSmith 是的，但关键是得到正确的结果，模数 2 ** 48。这就是 uint32 的工作原理，例如：在 C 中将两个 uint32 相乘时也会丢失精度，但结果仍然是正确的模数 2 ** 32。错误的u48_mul 实际上会给你错误的结果，相差 1 或 2。

这很棘手，因为 Javascript 二元运算符会将其操作数转换为 32 位整数。

【参考方案1】：

假设 a, b >= 0，如果您将 a 和 b 写为基数 2²⁴ 整数，则您有 a = a₁⋅2²⁴ + a₀ 和 b = b₁⋅2^{24 + b₀, 0 1, a₀, b₁, b₀24支持>.}

在这种形式中 a⋅b = a₁⋅b₁⋅2⁴⁸ + (a₁ ⋅b₀ + a₀⋅b₁)⋅2²⁴ + a₀⋅b₀。现在取这个 mod 2⁴⁸ 会使第一项变为 0。像 a₁⋅b₀ 这样的乘积是两个 24 位的乘积整数，因此乘以 JS 数字会产生精确的 48 位结果。将两个这样的值相加可能会产生一个 49 位的值，但这仍然是 53 并且因此是准确的。由于我们要乘以 2²⁴，我们只需要保留这个和的低 24 位。这很好，因为当我们使用 24 位掩码进行 AND (&) 时，JS 将只保留低位 32 位。最后，我们将结果添加到另一个 48 位值 a₀⋅b₀。结果可能会超过 2⁴⁸，如果超过了，我们只需减去 2⁴⁸。

function mulmod48(a, b) 
    const two_to_the_24th = 16777216;
    const two_to_the_48th = two_to_the_24th * two_to_the_24th;
    const mask24 = two_to_the_24th - 1;

    let a0 = a & mask24;
    let a1 = (a - a0) / two_to_the_24th;
    let b0 = b & mask24;
    let b1 = (b - b0) / two_to_the_24th;
    let t1 = ((a1 * b0 + a0 * b1) & mask24) * two_to_the_24th;
    let result = t1 + a0 * b0;
    if (result >= two_to_the_48th) 
        result -= two_to_the_48th;
    
    return result;

感谢@Mark Dickinson 的观察，这可以通过将其中一个除法消除 2²⁴ 来改善。尽管 IEEE 754 binary64 浮点数可以准确表示 -2⁵³ 和 2⁵³ 之间的所有整数，但这些并不是唯一可以准确表示的整数。特别是在指定的 IEEE 754 指数范围内乘以 2 的幂的 48 或 49 位整数也可以精确表示。因此我们可以替换这五行

    let a0 = a & mask24;
    let a1 = (a - a0) / two_to_the_24th;
    let b0 = b & mask24;
    let b1 = (b - b0) / two_to_the_24th;
    let t1 = ((a1 * b0 + a0 * b1) & mask24) * two_to_the_24th;

这三个

    let a0 = a & mask24;
    let b0 = b & mask24;
    let t1 = ((((a - a0) * b0 + a0 * (b - b0)) / two_to_the_24th) & mask24) * two_to_the_24th;

因为 (a - a0) 和 (b - b0) 都是 2²⁴ 的倍数，所以它们与 b0 和 a0 的乘积也必须如此，因此这些乘积之和也必须如此。换句话说，(a - a0) 和 (b - b0) 都只有 24 个 significant 位，所以乘积只有 48 个有效位，它们的和只有 49 个有效位。

现在，这比使用 Bigints 更快吗？在Chrome 96.0.4664.55 (Official Build) (x86_64) 的实验中，它比您的 BigInt 版本的u48_mul 快了近 6 倍。

【讨论】：

哦，我只是想知道我是否缺少一些简单的解决方案。我现在为你在回答我的问题时付出的大量工作感到难过。这是非常令人印象深刻的。希望我能提示这个答案。谢谢。

我喜欢解决这些小问题，感谢您提出问题。

@MaiaVictor "希望我能给出这个答案" --> wish granted.

好吧，a1 * b0 + a0 * b1 将照常使用 IEEE 754 算法执行。它确实可能大于 2**48，但仍会精确计算，因为 a1 和 b1 并且结果都是 2**24 的精确整数倍。所以溢出不会导致错误的结果。

是的，乘积可以和2**72一样大，但它是2**24的整数倍，所以它只有48个有效位，并且是完全可以用 JavaScript 用于其Number 类型的 IEEE 754 binary64 格式表示。 IOW，所有算术仍然是精确的，最终结果将是完全正确的。