计算质数时的效率 Haskell

Question

我有以下一组函数来计算 Haskell 中小于或等于数字 n 的素数的数量。

该算法接受一个数字，检查它是否可以被 2 整除，然后检查它是否可以被奇数整除，直到被检查的数字的平方根为止。

-- is a numner, n, prime? 
isPrime :: Int -> Bool
isPrime n = n > 1 &&
              foldr (\d r -> d * d > n || (n `rem` d /= 0 && r))
                True divisors

-- list of divisors for which to test primality
divisors :: [Int]
divisors = 2:[3,5..]

-- pi(n) - the prime counting function, the number of prime numbers <= n
primesNo :: Int -> Int
primesNo 2 = 1
primesNo n
    | isPrime n = 1 + primesNo (n-1)
    | otherwise = 0 + primesNo (n-1)

main = print $ primesNo (2^22)

使用带有 -O2 优化标志的 GHC，在我的系统上计算 n = 2^22 的素数需要大约 3.8 秒。以下 C 代码耗时约 0.8 秒：

#include <stdio.h>
#include <math.h>

/*
    compile with: gcc -std=c11 -lm -O2 c_primes.c -o c_orig
*/

int isPrime(int n) {
    if (n < 2)
        return 0;
    else if (n == 2)
        return 1;
    else if (n % 2 == 0)
        return 0;
    int uL = sqrt(n);
    int i = 3;
    while (i <= uL) {
        if (n % i == 0)
            return 0;
        i+=2;
    }
    return 1;
}

int main() {
    int noPrimes = 0, limit = 4194304;
    for (int n = 0; n <= limit; n++) {
        if (isPrime(n))
            noPrimes++;
    }
    printf("Number of primes in the interval [0,%d]: %d\n", limit, noPrimes);
    return 0;
}

该算法在 Java 中大约需要 0.9 秒，在 JavaScript（在节点上）中大约需要 1.8 秒，所以感觉 Haskell 版本比我预期的要慢.无论如何，我可以在 Haskell 中更有效地编码而不改变算法吗？

编辑

@dfeuer 提供的以下版本的 isPrime 将运行时间缩短了一秒，将其降至 2.8 秒（从 3.8 秒下降）。尽管这仍然比 JavaScript (Node) 慢，后者需要大约 1.8 秒，如此处所示，Yet Another Language Speed Test.

isPrime :: Int -> Bool
isPrime n
  | n <= 2 = n == 2
  | otherwise = odd n && go 3
  where
    go factor
      | factor * factor > n = True
      | otherwise = n `rem` factor /= 0 && go (factor+2)

编辑

在上面的isPrime函数中，函数go调用了factor * factor 对于单个 n 的每个除数。我想将 factor 与 n 的平方根进行比较会更有效，因为每个 [=] 只需计算一次32=]n。但是，使用下面的代码，计算时间增加了大约 10%，是不是每次计算不等式时都会重新计算 n 的平方根（对于每个因子)?

isPrime :: Int -> Bool
isPrime n
  | n <= 2 = n == 2
  | otherwise = odd n && go 3
  where
    go factor
      | factor > upperLim = True
      | otherwise = n `rem` factor /= 0 && go (factor+2)
      where
        upperLim = (floor.sqrt.fromIntegral) n

Answer 1

我强烈建议您使用不同的算法，例如 paper by Melissa O'Neill, or the version used in Math.NumberTheory.Primes from the arithmoi 包中讨论的埃拉托色尼筛法，它还提供优化的素数计数功能。但是，这可能会让您获得更好的常数因子：

-- is a number, n, prime? 
isPrime :: Int -> Bool
isPrime n
  | n <= 2 = n == 2
  | otherwise = odd n && -- Put the 2 here instead
        foldr (\d r -> d * d > n || (n `rem` d /= 0 && r))
                True divisors

-- list of divisors for which to test primality
divisors :: [Int]
{-# INLINE divisors #-} -- No guarantee, but it might possibly inline and stay inlined,
               -- so the numbers will be generated on each call instead of
               -- being pulled in (expensively) from RAM.
divisors = [3,5..] -- No more 2:

摆脱 2: 的原因是称为 "foldr/build fusion"、"short cut deforestation" 或 "list fusion" 的优化可能会使您的除数列表消失离开，但是，至少在 GHC < 7.10.1 的情况下，2: 将阻止优化。

编辑：这似乎不适合你，所以这里有一些其他的尝试：

isPrime n
  | n <= 2 = n == 2
  | otherwise = odd n && go 3
  where
    go factor
      | factor * factor > n = True
      | otherwise = n `rem` factor /= 0 && go (factor+2)

Answer 2

内联一切，松开多余的测试，添加严格性注释以确保：

{-# LANGUAGE BangPatterns #-}

-- pi(n) - the prime counting function, the number of prime numbers <= n
primesNo :: Int -> Int
primesNo n
    | n < 2 = 0
    | otherwise = g 3 1
 where
   g  k !cnt | k > n     = cnt
             | go 3      = g  (k+2) (cnt+1)
             | otherwise = g  (k+2) cnt
      where go f 
               | f*f > k   = True
               | otherwise = k `rem` f /= 0 && go (f+2) 

main = print $ primesNo (2^22)

go 测试功能与 dfeuer 的回答相同。像往常一样使用 -O2 进行编译，并且 始终通过运行独立的可执行文件 （使用类似 > test +RTS -s 的东西）进行测试。

可以直接调用 g（这确实是微优化）：

primesNo n
    | n < 2 = 0
    | otherwise = g 3 1
 where
   g  k !cnt | k > n     = cnt
             | otherwise = go 3
      where go f 
               | f*f > k        = g (k+2) (cnt+1)
               | k `rem` f == 0 = g (k+2) cnt
               | otherwise      = go (f+2)

可能会或可能不会加快速度的更实质性的变化（仍然保持算法可以说是相同的）是将其翻转过来，以节省平方计算：通过 [3] 从 9 开始的所有赔率进行测试到 23，通过 [3,5] 从 25 到 47 的所有赔率，等等，沿着 this segmented code:

import Data.List (inits)

primesNo n = length (takeWhile (<= n) $ 2 : oddprimes)
  where
    oddprimes = sieve 3 9 [3,5..] (inits [3,5..]) 
    sieve x q ~(_:t) (fs:ft) =
      filter ((`all` fs) . ((/=0).) . rem) [x,x+2..q-2]
      ++ sieve (q+2) (head t^2) t ft

有时将代码调整为使用 and 而不是 all 也会改变速度。可以通过内联和简化一切来尝试进一步加速（用计数等替换 length）。

Answer 3

一般来说，我发现 Haskell 中的循环比用 C 可以完成的循环要慢 3-4 倍。

为了帮助理解性能差异，我稍微修改了程序，以便每次迭代进行固定数量的除数测试并添加了一个参数 e 来控制进行了多少次迭代 - 执行的（外部）迭代次数为 2^e。对于每次外部迭代约进行了 2^21 个除数测试。

每个程序和脚本的源代码运行并分析结果可在此处找到：https://github.com/erantapaa/loopbench

欢迎提出改进基准测试的请求。

这是我使用 ghc 7.8.3（在 OSX 下）在 2.4 GHz Intel Core 2 Duo 上获得的结果。使用的 gcc 是 "Apple LLVM version 6.0 (clang-600.0.56) (based on LLVM 3.5svn)".

e     ctime   htime  allocated  gc-bytes alloc/iter  h/c      dns
10   0.0101  0.0200      87424      3408             1.980   4.61
11   0.0151  0.0345     112000      3408             2.285   4.51
12   0.0263  0.0700     161152      3408             2.661   5.09
13   0.0472  0.1345     259456      3408             2.850   5.08
14   0.0819  0.2709     456200      3408             3.308   5.50
15   0.1575  0.5382     849416      9616             3.417   5.54
16   0.3112  1.0900    1635848     15960             3.503   5.66
17   0.6105  2.1682    3208848     15984             3.552   5.66
18   1.2167  4.3536    6354576     16032  24.24      3.578   5.70
19   2.4092  8.7336   12646032     16128  24.12      3.625   5.75
20   4.8332 17.4109   25229080     16320  24.06      3.602   5.72

e          = exponent parameter
ctime      = running time of the C program
htime      = running time of the Haskell program
allocated  = bytes allocated in the heap (Haskell program)
gc-bytes   = bytes copied during GC (Haskell program)
alloc/iter = bytes allocated in the heap / 2^e
h / c      = htime divided by ctime
dns        = (htime - ctime) divided by the number of divisor tests made
             in nanoseconds

# divisor tests made = 2^e * 2^11

一些观察：

Haskell 程序以每次（外）循环迭代大约 24 字节的速率执行堆分配。 C 程序显然不执行任何分配，运行s 完全在 L1 缓存中。
对于 10 到 14 之间的 e，gc-bytes 计数保持不变，因为没有对那些运行s 执行垃圾收集。
时间比率 h/c 随着分配的增加而逐渐变差。
dps 是 Haskell 程序每次除数测试所花费的额外时间的度量；它随着分配总量的增加而增加。还有一些高原表明这是由于缓存效应。

众所周知，GHC 不会生成与 C编译器产生。您支付的罚款约为。每次迭代 4.6 ns。此外，看起来 Haskell 也受到缓存效应的影响，因为堆分配。

每次分配 24 字节和每次循环迭代 5 ns 对于大多数程序，但是当您有 2^20 次分配和 2^40 次循环迭代时它成为一个因素。

Answer 4

C代码使用32位整数，而Haskell代码使用64位整数。

原始 C 代码在我的计算机上运行 0.63 秒。但是，如果我用 long-s 替换 int-s，它会在 2.07 秒内运行 gcc 和 2.17 秒与 clang.

相比之下，更新后的 isPrime 函数（在线程问题中查看）运行时间为 2.09 秒（使用 -O2 和 -fllvm）。请注意，这比 clang 编译的 C 代码略好，即使它们使用相同的 LLVM 代码生成器。

原始Haskell代码运行3.2秒，我认为为了方便使用列表进行迭代，这是可以接受的开销。

计算质数时的效率 Haskell

Efficiency in Haskell when counting primes

performance

primes

haskell

编辑

编辑