23

我在考虑一个大数除法的算法：用余数除以 bigint C 除以 bigint D，我们知道 C 在基数 b 中的表示，而 D 的形式为 b^k-1。在示例中显示它可能是最容易的。让我们尝试将 C=21979182173 除以 D=999。

  • 我们将数字写成三位数：21 979 182 173
  • 我们取连续集合的总和（模 999），从左边开始：21 001 183 356
  • 我们将 1 添加到我们“超过 999”之前的那些集合：22 001 183 356

实际上，21979182173/999=22001183 和余数 356。

我已经计算了复杂度，如果我没记错的话，算法应该在 O(n) 中工作，n 是基 b 表示中 C 的位数。我还在 C++ 中做了一个非常粗略和未优化的算法版本（仅适用于 b=10），针对 GMP 的通用整数除法算法对其进行了测试，它确实似乎比 GMP 更好。我在任何地方都找不到这样的实现，所以我不得不求助于对一般部门进行测试。

我发现几篇文章讨论了似乎非常相似的问题，但没有一篇文章专注于实际实现，特别是在不同于 2 的基础上。我想这是因为数字在内部存储的方式，尽管提到的算法似乎有用，比如说，b = 10，即使考虑到这一点。我也尝试联系其他人，但同样无济于事。

因此，我的问题是：是否有文章或书籍或其他东西描述了上述算法，可能讨论了实现？如果不是，那么我尝试在 C/C++ 中实现和测试这样的算法是否有意义，或者这个算法在某种程度上是天生不好的？

另外，我不是程序员，虽然我在编程方面相当不错，但我承认我对计算机“内部结构”知之甚少。因此，请原谅我的无知——这篇文章中很可能有一件或多件非常愚蠢的事情。再次抱歉。

非常感谢！

进一步澄清评论/答案中提出的观点：

谢谢大家——因为我不想用同样的东西评论所有很棒的答案和建议，我只想谈谈你们中很多人提到的一点。

我完全清楚，一般来说，在 2^n 基础上工作显然是最有效的做事方式。几乎所有 bigint 库都使用 2^32 或其他。但是，如果（而且，我强调，它只对这个特定的算法有用！）我们将 bigints 实现为以 b 为底的数字数组，该怎么办？当然，我们在这里要求 b 是“合理的”：b=10，最自然的情况，似乎足够合理。我知道考虑到内存和时间，考虑到数字是如何在内部存储的，我知道它或多或少是低效的，但我已经能够，如果我的（基本的和可能有某种缺陷的）测试是正确的，比 GMP 的一般部门更快地产生结果，这将对实现这样的算法有意义。

Ninefingers 注意到在这种情况下我必须使用昂贵的模运算。我希望不会：我可以通过查看 old+new+1 的位数来查看 old+new 是否交叉，例如 999。如果它有 4 位数字，我们就完成了。更重要的是，由于old<999和new<=999，我们知道如果old+new+1有4位（不能多），那么，(old+new)%999等于删除( old+new+1)，我认为我们可以便宜地做到这一点。

当然，我不是在争论这个算法的明显局限性，也不是我声称它不能被改进——它只能除以特定类别的数字，我们必须先验地知道以 b 为基数的被除数的表示。但是，例如，对于 b=10，后者似乎很自然。

现在，假设我们已经实现了我上面概述的 bignums。说以 b 为底的 C=(a_1a_2...a_n) 和 D=b^k-1。该算法（可能会更加优化）会像这样。希望没有太多错别字。

  • 如果k>n，我们显然完成了
  • 在 C 的开头添加一个零（即 a_0=0）（以防我们尝试将 9999 与 99 相除）
  • l=n%k （“常规”整数的 mod - 不应该太贵）
  • old=(a_0...a_l) （第一组数字，可能少于 k 位）
  • for (i=l+1; i < n; i=i+k) （我们将有 floor(n/k) 次左右的迭代）
    • 新=(a_i...a_(i+k-1))
    • new=new+old （这是 bigint 加法，因此 O(k)）
    • aux=new+1 （再次，bigint 加法 - O(k) - 我不满意）
    • 如果 aux 有超过 k 个数字
      • 删除 aux 的第一个数字
      • old=old+1 (bigint 再次加法)
      • 在开头用零填充旧的，所以它应该有尽可能多的数字
      • (a_(ik)...a_(i-1))=旧（如果 i = l+1，( a_0 ... a_l )=旧）
      • 新=辅助
    • 在开头用零填充新的，所以它应该有尽可能多的数字
    • (a_i...a_(i+k-1)=新
  • quot=(a_0...a_(n-k+1))
  • rem=新

在那里，感谢您与我讨论这个问题 - 正如我所说，这在我看来确实是一个有趣的“特例”算法，如果没有人看到它有任何致命缺陷，可以尝试实施、测试和讨论。如果这是迄今为止尚未广泛讨论的事情，那就更好了。请让我知道你的想法。对不起，很长的帖子。

此外，还有一些个人评论：

@Ninefingers：我实际上对 GMP 的工作原理、它的作用以及一般的 bigint 除法算法有一些（非常基本的！）知识，所以我能够理解你的大部分论点。我也知道 GMP 是高度优化的，并且在某种程度上为不同的平台定制了自己，所以我当然不会试图“打败它”——这似乎与用尖头棍子攻击坦克一样富有成效。然而，这不是这个算法的想法——它适用于非常特殊的情况（GMP 似乎没有涵盖）。在不相关的说明中，您确定在 O(n) 中完成了一般划分吗？我见过的最多的是M（n）。（如果我理解正确，那在实践中（Schönhage-Strassen 等）可能不会达到 O(n)。如果我是正确的，Fürer 的算法仍然没有达到 O(n)，几乎纯粹是理论上的。）

@Avi Berger：尽管想法相似，但这实际上似乎与“淘汰九点”并不完全相同。但是，如果我没记错的话，上述算法应该一直有效。

4

3 回答 3

12

您的算法是一种以 10 为基数的算法的变体，称为“投出 9”。您的示例使用基数 1000 和“淘汰”999（比基数少一个）。这曾经在小学教过，作为快速检查手算的方法。我有一位高中数学老师，得知它不再被教了，吓坏了，并填补了我们的空白。

在基数 1000 中剔除 999 不能作为通用除法算法。它将生成与实际商和余数模 999 一致的值 - 而不是实际值。你的算法有点不同，我没有检查它是否有效，但它基于有效地使用基数 1000 并且除数比基数小 1。如果您想尝试除以 47，则必须先转换为以 48 为基数的数字系统。

谷歌“投出九分”以获取更多信息。

编辑：我最初读你的帖子有点太快了，你确实知道这是一种有效的算法。正如@Ninefingers 和@Karl Bielefeldt 在他们的评论中比我更清楚地陈述的那样，您在绩效评估中没有包括的是转换为适合手头特定除数的基数。

于 2011-02-23T22:54:07.403 回答
5

我觉得有必要根据我的评论对此进行补充。这不是答案，而是对背景的解释。

bignum 库使用所谓的肢体——在 gmp 源中搜索 mp_limb_t，它通常是一个固定大小的整数字段。

当您执行加法之类的操作时，一种方法（尽管效率低下）是这样做：

doublelimb r = limb_a + limb_b + carryfrompreviousiteration

在总和大于肢体大小的情况下，这个双倍大小的肢体捕获了肢体_a + 肢体_b 的溢出。因此，如果总数大于 2^32，如果我们使用 uint32_t 作为肢体大小，则可以捕获溢出。

我们为什么需要这个？好吧，你通常做的是循环遍历所有的肢体——你已经自己完成了将整数除以并遍历每个肢体——但我们首先做 LSL（所以首先是最小的肢体），就像你做算术一样用手。

这可能看起来效率低下，但这只是 C 的做事方式。要真正突破大炮，x86 有adc一个指令 - 加进位。这是一个算术和你的字段，如果算术溢出寄存器的大小，则设置进位位。下次执行addoradc时，处理器也会考虑进位位。在减法中，它被称为借用标志。

这也适用于班次操作。因此，处理器的这一特性对于使 bignums 快速运行至关重要。所以事实是，芯片中有电子电路来做这些事情——用软件做总是会更慢。

无需过多介绍，操作就是通过这种加、移、减等能力建立起来的。它们是至关重要的。哦，如果你做得对，你会使用每个肢体的处理器寄存器的整个宽度。

第二点——碱基之间的转换。您不能在数字中间取值并更改其基数，因为您无法解释原始基数中其下方数字的溢出，并且该数字无法解释下方数字的溢出。 .. 等等。简而言之，每次你想改变 base 时，你都需要将整个 bignum从原来的 base 重新转换回你的新 base。所以你必须至少走三遍（所有的四肢）。或者，在所有其他操作中昂贵地检测溢出......记住，现在你需要做模运算来计算你是否溢出，而在处理器为我们做这件事之前。

我还想补充一点，虽然在这种情况下你所拥有的可能很快，但请记住，作为一个 bignum 库，gmp 为你做了很多工作，比如内存管理。如果您正在使用mpz_，那么您正在使用我在此处描述的抽象，对于初学者。最后，gmp 对您听说过的几乎所有平台以及更多平台都使用了带有展开循环的手动优化组装。它与 Mathematica、Maple 等人一起提供是有充分理由的。

现在，仅供参考，一些阅读材料。

  • 现代计算机算术是针对任意精度库的类似 Knuth 的工作。
  • Donald Knuth，半数值算法（计算机编程艺术第二卷）。
  • William Hart关于实现bsdnt算法的博客，他在其中讨论了各种除法算法。如果您对 bignum 库感兴趣，这是一个很好的资源。在我开始关注这类东西之前，我认为自己是一个优秀的程序员......

总结一下：除法汇编指令很烂，所以人们通常计算逆和乘法，就像你在模算术中定义除法时所做的那样。现有的各种技术（参见 MCA）大多是 O(n)。

编辑：好的，并不是所有的技术都是 O(n)。大多数称为 div1 的技术（除以不大于肢体的东西是 O(n)。当你变大时，你最终会得到 O(n^2) 复杂度；这很难避免。

现在，您可以将 bigints 实现为数字数组吗？嗯，是的，你当然可以。但是，考虑一下加法下的想法

/* you wouldn't do this just before add, it's just to 
   show you the declaration.
 */
uint32_t* x = malloc(num_limbs*sizeof(uint32_t));
uint32_t* y = malloc(num_limbs*sizeof(uint32_t));
uint32_t* a = malloc(num_limbs*sizeof(uint32_t));
uint32_t m;

for ( i = 0; i < num_limbs; i++ )
{
    m = 0;
    uint64_t t = x[i] + y[i] + m;
    /* now we need to work out if that overflowed at all */
    if ( (t/somebase) >= 1 ) /* expensive division */
    {
        m = t % somebase; /* get the overflow */
    }
}

/* frees somewhere */

这是您希望通过您的方案添加的内容的粗略草图。所以你必须在碱基之间进行转换。因此，您将需要转换为基础的表示形式，然后在完成后返回，因为这种形式在其他任何地方都非常慢。我们在这里不是在讨论 O(n) 和 O(n^2) 之间的区别，而是在讨论每个肢体的昂贵除法指令或每次您想要除法时的昂贵转换。看到这个。

接下来，您如何将您的部门扩展到一般案件部门？我的意思是当你想从上面的代码中除以这两个数字x和y时。答案是，如果不求助于昂贵的基于 bignum 的设施，你就无法做到。见克努特。取模数大于您的大小是行不通的。

让我解释。试试 21979182173 mod 1099。为了简单起见，我们在这里假设我们可以拥有的最大尺寸字段是三位数。这是一个人为的例子，但我知道的最大字段大小是否使用 128 位使用 gcc 扩展。无论如何，重点是，你：

21 979 182 173

把你的号码分成四肢。然后你取模和求和：

21 1000 1182 1355

它不起作用。这是 Avi 正确的地方，因为这是一种排除 9 的形式，或者它的改编形式，但它在这里不起作用，因为我们的字段一开始就溢出 - 您正在使用模数来确保每个字段都保持在范围内它的肢体/字段大小。

那么解决方案是什么？把你的号码分成一系列大小合适的大数字？并开始使用 bignum 函数来计算您需要的一切？这将比任何现有的直接操作字段的方式慢得多。

现在，也许您只是提出这种情况来除以肢体，而不是大数，在这种情况下它可以工作，但是亨塞尔除法和预先计算的逆等在没有转换要求的情况下可以做到。我不知道这个算法是否会比亨塞尔除法更快；这将是一个有趣的比较；问题来自bignum 库中的通用表示。在现有的 bignum 库中选择的表示是出于我已经扩展的原因 - 它在汇编级别是有意义的，它是第一次完成的地方。

作为旁注；你不必用它uint32_t来代表你的四肢。您最好使用系统寄存器的大小（例如 uint64_t），以便您可以利用汇编优化版本。因此，在 64 位系统上adc rax, rbx，如果结果溢出 2^64 位，则仅设置溢出 (CF)。

tl；博士版本：问题不在于您的算法或想法；这是在基数之间转换的问题，因为您的算法所需的表示并不是在 add/sub/mul 等中最有效的方法。套用 knuth 的话说：这向您展示了数学优雅和计算效率之间的区别。

于 2011-02-24T00:05:45.330 回答
0

如果您需要经常除以相同的除数，则使用它（或它的幂）作为您的基数使得除法与位移对于基数 2 二进制整数一样便宜。

如果你愿意，你可以使用基数 999；使用 10 的幂并没有什么特别之处，只是它使转换为十进制整数非常便宜。（您可以一次工作一个肢体，而不必对整个整数进行完全除法。这就像将二进制整数转换为十进制与将每 4 位转换为十六进制数字之间的区别。二进制 -> 十六进制可以开始具有最高有效位，但转换为非 2 次幂的基必须是 LSB 优先使用除法。）

例如，要计算具有性能要求的代码高尔夫问题的 Fibonacci(10 9 ) 的前 1000 个十进制数字，我的 105 字节 x86 机器代码答案使用了与此 Python 答案相同的算法：通常的a+=b; b+=aFibonacci 迭代，但是每次除以（的幂）10a变得太大。

斐波那契的增长速度比进位传播的速度更快，因此偶尔丢弃低位小数位并不会长期改变高位位。（你保留一些超出你想要的精度的额外部分）。

除以 2 的幂是行不通的，除非您跟踪您丢弃了多少 2 的幂，因为最后的二进制 -> 十进制转换将取决于此。

所以对于这个算法，你必须做扩展精度加法，然后除以 10（或者你想要的任何 10 的幂）。

我将基数为 10 的9 个肢体存储在 32 位整数元素中。除以 10 9非常便宜：只是一个指针增量来跳过下肢。我没有实际执行 a memmove，而是偏移了下一次添加迭代使用的指针。

我认为除以 10 ^ 9 以外的 10 的幂会有点便宜，但需要在每个肢体上进行实际除法，并将余数传播到下一个肢体。

这种方式的扩展精度加法比使用二进制肢体要贵一些，因为我必须使用比较手动生成进位：（sum[i] = a[i] + b[i]; carry = sum < a;无符号比较）。并且还使用条件移动指令根据该比较手动换行到 10^9。但是我能够使用该进位作为adc（x86 add-with-carry 指令）的输入。

你不需要一个完整的模来处理加法的包装，因为你知道你最多包装一次。

这浪费了每个 32 位肢体的 2 位多一点：10^9 而不是2^32 = 4.29... * 10^9. 每个字节存储一个以 10 为基数的数字会显着降低空间效率，并且对性能非常不利，因为 8 位二进制加法的成本与现代 64 位 CPU 上的 64 位二进制加法相同。

我的目标是代码大小：为了纯粹的性能，我会使用持有 base-10^19“数字”的 64 位肢体。( ，因此每 64 位浪费了不到 1 位。) 这使您可以使用每条硬件指令2^64 = 1.84... * 10^19完成两倍的工作。add嗯，实际上这可能是一个问题：两个分支的总和可能包含 64 位整数，因此仅检查> 10^19是不够的。您可以在 base5*10^18或 base 中工作10^18，或者进行更复杂的进位检测，检查二进制进位和手动进位。

以每 4 位半字节一个数字存储压缩 BCD 对性能来说会更糟，因为没有硬件支持阻止一个字节内从一个半字节到下一个半字节的进位。

总体而言，我的版本在相同硬件上的运行速度比 Python 扩展精度版本快大约 10 倍（但它有很大的速度优化空间，通过减少划分频率）。（70 秒或 80 秒对比 12 分钟）

不过，我认为对于该算法的这种特定实现（我只需要加法和除法，并且每隔几次加法就会发生除法），选择 base-10^9 肢体非常好。对于第 N 个斐波那契数，有更高效的算法，不需要进行 10 亿次扩展精度加法。

于 2017-11-16T21:41:11.907 回答