您的Haskell仅（非）阶乘且适合

当谈到喜欢的语言时，我通常会说，在所有条件都相同的情况下，我更喜欢C ++作为数字粉碎机，而Haskell则用于其他所有事物。 定期检查该除法是否合理是有用的，直到最近才出现一个闲置且非常简单的问题：一个数字的所有除数之和如何随该数字的增长而变化（例如，对于前十亿个数字）。 威吓这项任务很容易（称其为结果磨床者是可耻的），因此，对于这种检查来说，这似乎是一个不错的选择。

此外，我仍然无法准确预测Haskell代码的性能，因此尝试尝试一些已知的不良方法以查看性能如何下降非常有用。

好吧，此外，您可以轻松地展示一种比对每个数字的除数进行正面搜索的更有效的算法 $$$1$ 之前 $$$n$ 。

演算法

因此，让我们从算法开始。

如何找到一个数字的所有除数之和 $$$n$ ？ 你可以经历所有 $$$k_1 \ in \ {1 \点\ lflo \ sqrt n \ rfloor \}$ 而对于每一个这样的 $$$k_1$ 检查除法的其余部分 $$$n$ 在 $$$k_1$ 。 如果剩下的是 $$$0$ ，然后添加到电池 $$$k_1 + k_2$ 在哪里 $$$k_2 = \ frac {n} {k_1}$ 如果 $$$k_1 \ neq k_2$ ，而 $$$k_1$ 否则。

该算法可以应用吗 $$$n$ 时间，对于每个数字 $$$1$ 之前 $$$n$ ？ 当然可以。 会有什么困难？ 容易看到那个顺序 $$$O（n ^ \ frac {3} {2}）$ 除法-对于每个数字，我们都精确地得出除法的根，并且我们有数字 $$$n$ 。 我们可以做得更好吗？ 事实证明是的。

这种方法的问题之一是我们浪费了太多的精力。 太多的分裂并不能使我们成功，剩下的非零。 尝试变得更懒惰并从另一侧着手完成任务是很自然的：让我们为除数生成各种候选，然后看看它们满足什么数字？

因此，现在让我们需要为每个 $$$1$ 之前 $$$n$ 计算所有除数的总和。 为此，请遍历所有 $$$k_1 \ in \ {1 \点\ lflo \ sqrt n \ rfloor \}$ ，并且每个 $$$k_1$ 让我们经历所有 $$$k_2 \ in \ {k_1 \点\ lfloor \ frac {n} {k} \ rfloor \}$ 。 每对 $$$（k_1，k_2）$ 添加到具有索引的单元格 $$$k_1 \ cdot k_2$ 价值 $$$k_1 + k_2$ 如果 $$$k_1 \ neq k_2$ 和 $$$k_1$ 否则。

这个算法确实做到了 $$$n ^ \ frac {1} {2}$ 除法，并且每次乘法（比除法便宜）都会使我们获得成功：在每次迭代中，我们都会增加一些东西。 这比正面方法更为有效。

此外，采用相同的正面方法，您可以比较这两种实现，并确保它们以相当小的数量给出相同的结果，这应该会增加一点信心。

第一次实施

顺便说一下，这几乎完全是Haskell初始实现的伪代码：

module Divisors.Multi(divisorSums) where import Data.IntMap.Strict as IM divisorSums :: Int -> Int divisorSums n = IM.fromListWith (+) premap IM.! n where premap = [ (k1 * k2, if k1 /= k2 then k1 + k2 else k1) | k1 <- [ 1 .. floor $ sqrt $ fromIntegral n ] , k2 <- [ k1 .. n `quot` k1 ] ] 

Main模块很简单，我不带。

另外，这里我们只显示最多的数量 $$$n$ 为了便于与其他实现进行比较。 尽管Haskell是一种惰性语言，但在这种情况下，所有金额都将被计算（尽管其充分依据超出了本文的范围），因此无法算出我们不会无意中计算任何数字。

它工作多快？ 在我的i7 3930k上，一个流在0.4 s内完成了100,000个元素的计算。 在这种情况下，在计算上花费0.15秒，在GC上花费0.25秒。 而且，由于int的大小为8字节，因此我们占用了大约8兆字节的内存，理想情况下，我们应该有800千字节。

好（不是真的）。 这些数字将如何随着数量的增加而增长？ 对于1 000 000个元素，它已经工作了大约7.5秒，在计算上花费了3秒，在GC上花费了4.5秒，并且还占用了80兆字节（超出所需数量的10倍）。 并且即使我们假装一秒钟成为高级Java软件开发人员并开始调整GC，我们也不会有太大改变。 不好 看来我们永远不会等待十亿个数字，也不会进入内存：我的计算机上只有64 GB的RAM，如果这种趋势继续下去，大约需要80 GB。

似乎该做点时间了

C ++选项

让我们尝试去理解什么是有意义的，为此，我们将在代码上编写代码。

好吧，由于我们已经有了调试的算法，所以一切都很简单：

 #include <vector> #include <string> #include <cmath> #include <iostream> int main(int argc, char **argv) { if (argc != 2) { std::cerr << "Usage: " << argv[0] << " maxN" << std::endl; return 1; } int64_t n = std::stoi(argv[1]); std::vector<int64_t> arr; arr.resize(n + 1); for (int64_t k1 = 1; k1 <= static_cast<int64_t>(std::sqrt(n)); ++k1) { for (int64_t k2 = k1; k2 <= n / k1; ++k2) { auto val = k1 != k2 ? k1 + k2 : k1; arr[k1 * k2] += val; } } std::cout << arr.back() << std::endl; } 

-O3 -march=native ，clang 8，一百万个元素在0.024 s内被处理，占用了分配的8 MB内存。 十亿-155秒，预期的8 GB内存。 哎呀 Haskell不好。 Haskell需要扔掉。 只有阶乘和原形就可以写了！ 还是不行

第二选择

显然，要通过IntMap运行所有生成的数据，也就是说，实际上是一个相对普通的地图-温和地说，而不是最明智的决定（是的，这与开始时提到的显然是糟糕的选择）。 为什么我们不使用C ++代码中的数组？

让我们尝试：

 module Divisors.Multi(divisorSums) where import qualified Data.Array.IArray as A import qualified Data.Array.Unboxed as A divisorSums :: Int -> Int divisorSums n = arr A.! n where arr = A.accumArray (+) 0 (1, n) premap :: A.UArray Int Int premap = [ (k1 * k2, if k1 /= k2 then k1 + k2 else k1) | k1 <- [ 1 .. floor bound ] , k2 <- [ k1 .. n `quot` k1 ] ] bound = sqrt $ fromIntegral n :: Double 

由于Int非常简单，在这里我们立即使用数组的未装箱版本，并且不需要惰性。 盒装版本仅在arr类型上有所不同，因此我们也不会迷失成语。 另外，这里bound是分开进行的，但这不是因为编译器是愚蠢的并且不执行LICM，而是因为您可以显式指定其类型，并避免编译器发出有关floor参数默认值的警告。

一百万个元素0.045秒（仅比加号差两倍！）。 8 MB的内存，在GC（！）中为零毫秒。 在较大的大小下，趋势仍然存在-比C ++慢大约两倍，并且内存量相同。 好结果！ 但是我们可以做得更好吗？

事实证明是的。 accumArray检查索引，在这种情况下我们不需要这样做-索引在构造上是正确的。 让我们尝试用accumArray替换对accumArray的调用：

 module Divisors.Multi(divisorSums) where import qualified Data.Array.Base as A import qualified Data.Array.IArray as A import qualified Data.Array.Unboxed as A divisorSums :: Int -> Int divisorSums n = arr A.! (n - 1) where arr = A.unsafeAccumArray (+) 0 (0, n - 1) premap :: A.UArray Int Int premap = [ (k1 * k2 - 1, if k1 /= k2 then k1 + k2 else k1) | k1 <- [ 1 .. floor bound ] , k2 <- [ k1 .. n `quot` k1 ] ] bound = sqrt $ fromIntegral n :: Double 

如您所见，除了需要从头开始编制索引（我认为这是库API中的错误，这是另一个问题）之外，所做的更改很小。 表现如何？

一百万个元素-0.021秒（哇，在误差范围内，但比专业人士还快！）。 自然，相同的8 MB内存，GC中相同的0 ms。

十亿个元素-152 s（看起来确实比加号还快！）。 略少于8 GB。 在GC中为0毫秒。 该代码仍然是惯用语。 我想我们可以说这是一场胜利。

总结

首先，令我感到惊讶的是，用unsafe版本替换accumArray会带来如此大的增长。 期望有10％到20％的比率是更合理的（毕竟，总的来说，用at()代替operator[]不会显着降低性能），但不能减半！

其次，在我看来，相当惯用的干净代码而不会产生可变的黏附达到这种性能水平，这真的很酷。

第三，当然，在所有级别上都可以进行进一步的优化。 举例来说，我敢肯定，您可以从附加代码中挤出更多。 但是，我认为，在所有此类基准测试中，在花费的工作量（和代码量）和最终的消耗之间取得平衡是很重要的。 否则，所有事情最终都会收敛到LLVM JIT或类似的挑战。 此外，肯定有更有效的算法可以解决此问题，但是此处提出的简短想法的结果也将适用于本周日的小冒险。

第四，我的最爱：类型系统需要开发。 这里不需要unsafe ，作为程序员，我可以证明循环中遇到的所有k_1, k_2 k_1 * k_2 <= n 。 在依赖类型语言的理想世界中，我将静态构造此证明并将其传递给相应的函数，这将消除运行时检查的需要。 但是，可惜，在Haskell中没有成熟的zavtipov，在有zavtipy（我知道）的语言中，也没有array和类似物。

第五，我对其他编程语言的了解还不足以满足这些语言的近基准测试要求，但是我的一位朋友用python写了一个类似物。 速度几乎快了一百倍，而内存则更糟。 而且算法本身非常简单，因此，如果有知识的人在注释中用Go，Rust，Julia，D，Java，Malbolge写一个类似物，或者共享一个比较，例如与C ++-他们机器上的代码-可能会很棒。

附言：对不起，请点击标题。 我无法提出更好的建议。