在 Haskell 中生成有限的素数列表

Question

在 Haskell 中有很多关于生成素数的主题，但在我看来，它们都依赖于 'isPrime' 函数，如果我们还不知道素数序列，应该如下所示：

isPrime k = if k > 1 then null [ x | x <- [2,3..(div k 2) + 1], k `mod` x == 0]
                     else False

（div 可能会替换为 sqrt，但仍然...）

我尝试根据“归纳定义”构造素数（假设我们有一组前 n 个素数，然后是 (n+1)th 个素数是最小整数，使得前 n 个素数都不是它的除数）。我试过用斐波那契数列的方式来做，即：

fib :: Int -> Int
fib 0 = 0
fib 1 = 1
fib n = fibs !! n
    where fibs = 0 : 1 : zipWith (+) fibs (tail fibs)

我最终得到了这个：

-- checking if second number is a divisor of first one
ifDoesn'tDivide :: Int -> Int -> Bool
ifDoesn'tDivide n k 
    | mod n k == 0 = False
    | otherwise    = True

-- generating list which consists of first n prime numbers
firstPrimes :: Int -> [Int]
-- firstPrimes 1  = [2]
firstPrimes n     = take n primes 
    where primes = 2:(tail primes) ++ 
         [head [x | x <- [3,4..], k <- primes, ifDoesn'tDivide x k == True]]

但它不起作用，堆栈溢出当n >= 2。关于如何修复它的任何建议？

“Haskell 可以根据自身定义数据结构，从而创建无限的数据结构”。前面提到的那些素数和斐波那契数列是根据自身定义数据结构的具体情况，斐波那契数列可以正常工作，但这些primes 不行。

我是否遗漏了什么，这两种算法在本质上是不同的吗？

附：所以，我想，我只是在寻找最“Haskellish”的方式来做到这一点。

Answer 1

你总是可以使用 Haskell 中相当优雅的筛子。

primes = sieve [2..]

sieve (p : xs) = p : sieve [ x | x <- xs, x `mod` p > 0 ]

所以要得到前 10 个素数

> take 10 primes
[2,3,5,7,11,13,17,19,23,29]

请注意，虽然isPrime 没有明确使用，但列表推导式确保列表中的每个数字都必须是相对于它之前的所有素数的素数，即素数。

这更有效，它是Eratosthenes' sieve的核心（编辑）。

上面的代码是第一个例子：

• 梅丽莎·E·奥尼尔，The Genuine Sieve of Eratosthenes

本文更详细地介绍了 Haskell 中筛子的有效实现以及惰性在计算中的作用。强烈推荐！

Answer 2

您的解决方案尝试的核心是primes 的无限列表，由：

primes = 2:(tail primes)
  ++ [head [x | x <- [3,4..]
              , k <- primes
              , ifDoesn'tDivide x k == True]]

更新：你在评论中提到你正在考虑这个算法，所以你想象 Haskell 将使用仍然为空的“当前”值tail primes，以便评估像[2] ++ [] ++ [3] 这样的东西，然后循环。但是，当然，Haskell 不是必须的，因此不能这样工作。在 Haskell 中，primes 有一个固定的定义，在整个程序执行过程中保持不变。 Haskell 程序可以逐渐“发现”（或更准确地说“计算”）定义，这允许我们首先根据自身定义primes，但它不能在执行过程中更改定义。

因此，在查看此定义时，您需要想象primes 和因此tail primes 在它们出现的任何地方都具有相同的值，即使在递归使用时也是如此。这与带参数的典型递归函数不同：

fact 0 = 1
fact n = n * fact (n-1)

这里，尽管 函数 fact 在它出现的任何地方都有相同的定义，左侧的 fact n 的值和右侧的 fact (n-1) 的值 -由于参数不同，hand side 可以不同。

无论如何，如果我们看这个primes 的定义，我们需要primes 是所有素数的无限列表它出现的所有地方（而不是一个改变或改变的值）随着时间的推移“增长”），那么你就会明白为什么这个定义不起作用了。在这里，primes 被定义为 2 : tail primes ++ [expr] 用于完成所有实际工作的复杂 expr，但 tail primes 应该是无限的，所以在评估这个表达式时，你甚至永远不会得到 em> 到expr，因为你永远不会用完列表tail primes。

即使忽略++ [expr] 位，因为primes 有一个固定定义，表达式如下：

primes = 2 : tail primes

不是根据自身定义无限列表的正确方法。问题是primes的第二个元素被定义为tail primes的第一个元素，也就是primes的第二个元素，所以primes的第二个元素被定义为自己。当 Haskell 试图“发现”/“计算”它的值时，这将创建一个无限循环。 fibs 定义的关键：

fibs = 0 : 1 : zipWith (+) fibs (tail fibs)

是第一个和第二个元素被给出，然后第三个元素被定义为前两个元素的总和——它不是根据自身定义的，而是根据前面的列表定义的元素。这是成功递归定义无限（甚至有限）列表的关键。

请注意，Haskell 并不“知道”fib 是一个无限列表，并且对无限列表没有任何特殊作用。这与递归定义的有限列表的工作方式相同：

countdown = 10 : takeWhile (> 0) (map (subtract 1) countdown)

关键还是countdown 中的每个元素的定义方式使其仅依赖于countdown 的previous 元素。

要修改您的primes 定义以这种方式工作，您可能想要做的是概括您​​的列表理解，从获取“2”之后的下一个素数到获得任何当前素数p 之后的“下一个”素数，基于primes 可用：

primeAfter p = head [x | x <- [p+1..], k <- primes, ifDoesn'tDivide x k]

这不起作用有两个原因。首先，因为primes 是无限的，这将永远检查不同k 值的可分性。我们需要修改它以检查素数 k 仅到当前素数 p：

primeAfter p = head [x | x <- [p+1..]
                       , k <- takeUntil (==p) primes
                       , ifDoesn'tDivide x k]

这使用了一个帮助器，它将列表的开头直到谓词为真的第一个元素：

takeUntil p lst = case break p lst of (a,y:b) -> a ++ [y]

其次，检查结构不正确。如果有 any 素数 k 不分割它，则此列表理解将允许通过 x。只有当 all 素数 k 不分割它时，我们才需要让 x 通过：

primeAfter p = head [x | x <- [p+1..]
                       , and [ifDoesn'tDivide x k
                               | k <- takeWhile (<=p) primes]]

那么它有机会工作，我们可以将primes定义为：

primes = go 2
  where go p = p : go (primeAfter p)

这里，go 将当前素数添加到列表中，然后使用primeAfter 递归到下一个素数。这是因为即使primeAfter p 访问由递归go 调用生成的无限列表primes，它也只使用该列表直到当前素数p，所以它只是停止在尝试访问列表中自己的值之前，仅使用在调用primeAfter p 之前生成的素数。

所以，这行得通，我认为这在很大程度上符合您最初尝试的精神：

-- note this simplified version:
ifDoesn'tDivide :: Int -> Int -> Bool
ifDoesn'tDivide n k = mod n k /= 0

primes :: [Int]
primes = go 2
  where go p = p : go (primeAfter p)

primeAfter :: Int -> Int
primeAfter p = head [x | x <- [p+1..]
                       , and [ifDoesn'tDivide x k
                               | k <- takeUntil (==p) primes]]

takeUntil :: (a -> Bool) -> [a] -> [a]
takeUntil p lst = case break p lst of (a,y:b) -> a ++ [y]

main :: IO ()
main = do
  print $ take 10 primes

正如@Mihalis 所指出的，primes 是 Haskell 中一个非常标准的示例，因此也有更优雅的单行解决方案。

Answer 3

TL;DR：不，这两种算法没有本质上的不同。

---

您的定义primes = 2:(tail primes) ++ .... 表示head primes = 2 和head (tail primes) = head ((tail primes) ++ ....) = head (tail primes)。这当然是有问题的，会导致无限递归。

---

在保留其意图的同时对代码进行最小修复可能是

firstPrimes1 :: Int -> [Int]
firstPrimes1 1  = [2]
firstPrimes1 n  = firstPrimes1 (n-1) ++ 
         take 1 [x | x <- [3,4..], 
                     and [ mod x k > 0 | k <- firstPrimes1 (n-1)]]

（这使用take 1 ... 代替您的[head ...]）。

速度慢得令人难以置信（looks 指数级，或更糟）。但它当然应该是，

firstPrimes2 1  = [2]
firstPrimes2 n  = let { ps = firstPrimes2 (n-1) } in
       ps ++ 
         take 1 [x | x <- [3,4..], 
                     and [ mod x k > 0 | k <- ps]]

现在非常慢，时间复杂度约为三次。但它应该真的是这样的，虽然：

firstPrimes2b 2  = [2]
firstPrimes2b n  = let { ps = firstPrimes2b (n-1) } in
       ps ++ 
         take 1 [x | x <- [last ps+1..], 
                     and [ mod x k > 0 | k <- ps]]

现在behaves 好像二次，而且在具体方面确实比它的前身快得多。

要像斐波那契流那样构造它，它可以写成

primes3 = 2 : concatMap foo [1..]
  where
  foo k = let { ps = take k primes3 } in
          take 1 [ x | x <- [last ps+1..], 
                       and [ mod x k > 0 | k <- ps]]
-- or 
primes4 = 2 : concatMap bar (tail (inits primes4))
  where
  bar ps = take 1 [ x | x <- [last ps+1..], 
                        and [ mod x k > 0 | k <- ps]]
-- or even 
primes5 = 2 : [p | (ps, q) <- zip (tail (inits primes5)) primes5
                 , p <- take 1 [ x | x <- [q+1..], 
                                     and [ mod x k > 0 | k <- ps]]]

确实，它看起来遵循一种归纳模式，特别是 complete 又名 "strong" 归纳，forall(n).(forall( k < n ).P(k)) => P(n)。

所以它与斐波那契计算没有根本不同，尽管后者仅指前两个元素，而这个指所有前面的元素，同时添加新的一个。但就像斐波那契流一样，这个序列最终也是根据自身定义的：primes = ..... primes ......。

inits 使bar 明确引用先前已知的素数ps，同时在每一步向它们添加一个（由take 1表示 )，就像你想要的那样。 concatMap 收集每次调用bar 产生的所有新的单元素段。

但为什么那应该只是 one 素数？我们不能从已知的先前素数k 安全地产生更多 个新素数吗？我们必须真的通过所有前面的素数来测试候选人，还是我们可以使用你在问题中也提到的众所周知的捷径？我们可以让它遵循完整的前缀归纳模式forall(n).(forall( k < floor(sqrt(n)) ).P(k)) => P(n)，这样只需要O(log log n)扩展步骤就可以到达第n个素数？

我们是否可以在每个步骤中从素数序列的每个前缀（当然，该序列始终保持相同）生成 更长 段，因此不指代每个候选者的所有前面的素数，但仅限于其中的一小部分？...

---

Eratosthenes 在 Haskell 中最直接表达的真正筛选是

import qualified Data.List.Ordered as O (minus)

primes = map head $ scanl (O.minus) [2..] [[p,p+p..] | p <- primes]

（minus 具有明显的语义，即使不从 data-ordlist 包中加载，也很容易自己实现。）

尽管 S. Horsley 牧师在 1772 年（重新？-）引入它时，^(*) 将 Eratosthenes 的筛子描述为相当于

oprimes = map head $ 
       scanl (O.minus . tail) [3,5..] [[p*p,p*p+2*p..] | p <- oprimes]

primes2 = 2 : oprimes

primesUpTo n = 2 : map head a ++ takeWhile (<= n) b
   where
   (a,b:_) = span ((<= n) . (^2) . head) $
       scanl (O.minus . tail) [3,5..] [[p*p,p*p+2*p..] | p <- oprimes]

运行length $ primesUpTo n 比length . takeWhile (<= n) primes 快得多。你知道为什么吗？

您能否修复primes2，使其在访问nth 元素时变得与primesUpTo 一样快？它可以按照您最初的想法，逐步扩展已知的素数段，如上一节所述。

另外，请注意这里根本没有使用isPrime 函数。这是埃拉托色尼的真筛子的标志，它不测试素数，它生成复合物，并免费获得复合物之间的素数。 p>

---

第一个scanl 代码如何工作：它以序列[2,3,4,5,...] 开头。然后它发出通知从其中删除[2,4,6,8,...]，并留下[3,5,7,9,...] 的等价物，即coprimes({2})。

^{（即使列表是无限的，这也是可行的，因为 Haskell 有惰性求值——只执行程序最终输出所需的计算量。）}

然后它发出通知从 them 列表中删除 [3,6,9,12,..]，并留下 coprimes({2,3})。

在每个阶段，它将head 从该时间点的序列中取出，并将该头元素放在一边，从而形成最终的素数序列。 p>

^{（同样可以用iterate（或unfoldr等）进行编码。这是一个很好的练习，可以帮助澄清那里到底发生了什么。当你这样做时，你会看到你将重新创建素数序列作为被迭代的 step 函数的参数的一部分（第一个 k 素数的共素数的当前序列，以及接下来，第 k+1 个素数，从该序列中删除 其 个倍数）。scanl 版本指的是 original 序列显式地取素数，一个接一个地从中取素数，但这是同一回事。）}

第二个scanl 变体仅枚举素数的奇数倍数，从素数的平方开始每个枚举（因此，例如3它是[9,15,21,27,...]，对于7 是[49,63,77,91,...])。它仍然开始枚举每个素数，而不是每个素数的平方。这就是为什么它必须在primesUpTo 函数中做出特殊安排，以便在可以停止时立即停止。这是its efficiency的关键。

---

^(*)pg 314 of Philosophical Transactions, Vol.XIII.

---

另请参阅：minus 已定义和使用 here，或 here。