为什么这个scala prime代这么慢/内存密集？

Question

在找到第10,001个素数时,我的内存耗尽.

object Euler0007 {
  def from(n: Int): Stream[Int] = n #:: from(n + 1)
  def sieve(s: Stream[Int]): Stream[Int] = s.head #:: sieve(s.filter(_ % s.head != 0))
  def primes = sieve(from(2))
  def main(args: Array[String]): Unit = {
    println(primes(10001))
  }
}

这是因为在每次"迭代"(这是这个上下文中的正确术语吗？)之后primes,我增加要调用的函数堆栈以使下一个元素一个接一个？

我在网上找到的一个解决方案是不使用迭代解决方案(我想避免进入函数式编程/惯用scala),这就是问题7(问题7):

lazy val ps: Stream[Int] = 2 #:: Stream.from(3).filter(i => ps.takeWhile(j => j * j <= i).forall(i % _ > 0))

从我所看到的,这不会导致这种类似递归的方式.这是一个很好的方法吗,或者你知道更好的方法吗？

Answer 1

这种缓慢的一个原因是它不是 Eratosthenes的筛子.阅读http://www.cs.hmc.edu/~oneill/papers/Sieve-JFP.pdf以获得详细说明(示例在Haskell中,但可以直接翻译成Scala).

我对Euler问题#7的旧解决方案也不是"真正的"筛子,但它似乎对于小数字而言足够好:

object Sieve {

    val primes = 2 #:: sieve(3)

    def sieve(n: Int) : Stream[Int] =
          if (primes.takeWhile(p => p*p <= n).exists(n % _ == 0)) sieve(n + 2)
          else n #:: sieve(n + 2)

    def main(args: Array[String]) {
      println(primes(10000)) //note that indexes are zero-based
    }
}

我认为你的第一个版本的问题是你只有defs和no val来收集结果,并且可以通过生成函数进行查询,所以你总是从头开始重新计算.

Answer 2

是的,这是因为你"增加要调用的函数堆栈以获取下一个元素,每次"迭代"之后一个 - 即每次获得每个素数后在过滤器堆栈顶部添加一个新过滤器.那个过滤器太多了.

这意味着每个生成的素数都会被其所有前面的素数进行测试 - 但只有那些低于其平方根的素数才真正需要.例如,为了获得10001个质数,104743,将有10000个创建过滤器,在运行时.但是下面只有66个素数323,平方根104743,所以实际上只需要66个滤波器.所有9934人都将在那里不必要地占据记忆,努力工作,完全没有附加价值.

这是"功能性筛子"的关键缺陷,它似乎起源于20世纪70年代的大卫特纳代码,后来又进入了SICP书和其他地方.它不是一个试验分区筛(而不是 Eratosthenes的筛子).这对它来说太遥远了.当最佳实施时,试验分割完全能够非常快地产生第10000个素数.

该代码的关键缺点是它不会在适当的时刻推迟过滤器的创建,最终会创建太多的过滤器.

现在谈论复杂性,"旧筛"代码是O(n ²),在n素数产生.最佳试验区分为O(n ^1.5/log ^0.5(n)),Eratosthenes筛为O(n*log(n)*log(log(n))).作为经验增长的顺序,第一个通常被视为~ n^2第二个~ n^1.45和第三个~ n^1.2.

您可以在此答案中找到基于Python生成器的代码,以实现最佳的试验分区(下半部分).它最初这里讨论对付哈斯克尔相当于你的筛功能.

就像一个例子,旧筛子的"可读伪代码" :)

primes = sieve [2..]  where
   sieve (x:xs) = x : sieve [ y | y <- xs, rem y x > 0 ]
                         -- list of 'y's, drawn from 'xs',
                         --         such that (y % x > 0)

并且对于最佳试验分区(TD)筛,在素数方块上同步,

primes = sieve [2..] primes   where
  sieve (x:xs) ps = x : (h ++ sieve [ y | y <- t, rem y p > 0 ] qs)
          where
            (p:qs) = ps     -- 'p' is head elt in 'ps', and 'qs' the rest
            (h,t)  = span (< p*p) xs    -- 'h' are elts below p^2 in 'xs'
                                        -- and 't' are the rest

和埃拉托色尼的筛子,由理查德·伯德设计,如在以另一种答案这里提到JFP文章,

primes = 2 : minus [3..] 
               (foldr (\p r-> p*p : union [p*p+p, p*p+2*p..] r) [] primes)
                      -- function of 'p' and 'r', that returns 
                      -- a list with p^2 as its head elt, ...

短而快 (minus a b是一个列表,a其中所有的elts b逐渐从中移除; union a b是一个列表,a其中所有elts b逐步添加到它中,没有重复;两者都处理有序的,非减少的列表).foldr是列表的正确折叠.由于它是线性的这个运行在~ n^1.33,以使其在运行~ n^1.2的树状折叠功能foldi都可以使用).

你的第二个问题的答案也是肯定的.你的第二个代码,用相同的"伪代码"重写,

ps = 2 : [i | i <- [3..], all ((> 0).rem i) (takeWhile ((<= i).(^2)) ps)]

与上面的最佳TD筛非常相似 - 两者都安排每个候选者通过其平方根以下的所有质数进行测试.虽然筛子使用推迟过滤器的运行时序列来排列,但后者的定义为每个候选者重新获取所需的素数.一个可能比另一个更快,具体取决于编译器,但两者基本相同.

第三个也是肯定的:Eratosthenes的筛子更好,

ps = 2 : 3 : minus [5,7..] (unionAll [[p*p, p*p+2*p..] | p <- drop 1 ps])

unionAll = foldi union' []          -- one possible implementation
union' (x:xs) ys = x : union xs ys
   -- unconditionally produce first elt of the 1st arg 
   -- to avoid run-away access to infinite lists

看起来它也可以在Scala中实现,从其他代码片段的相似性来判断.(虽然我不知道斯卡拉).unionAll这里实现了树状折叠结构(点击查看图片和完整代码),但也可以用滑动数组实现,沿着素数的多个流逐段工作.

TL; DR:是的,是的,是的.

Answer 3

FWIW,这是一个真正的Eratosthenes筛子:

def sieve(n: Int) = (2 to math.sqrt(n).toInt).foldLeft((2 to n).toSet) { (ps, x) => 
    if (ps(x)) ps -- (x * x to n by x) 
    else ps
}

这是一个无限的素数流,使用了Eratosthenes筛子的变体,保留了它的基本属性:

case class Cross(next: Int, incr: Int)

def adjustCrosses(crosses: List[Cross], current: Int) = {
  crosses map {
    case cross @ Cross(`current`, incr) => cross copy (next = current + incr)
    case unchangedCross                 => unchangedCross
  }
}

def notPrime(crosses: List[Cross], current: Int) = crosses exists (_.next == current)

def sieve(s: Stream[Int], crosses: List[Cross]): Stream[Int] = {
    val current #:: rest = s

    if (notPrime(crosses, current)) sieve(rest, adjustCrosses(crosses, current))
    else current #:: sieve(rest, Cross(current * current, current) :: crosses)
}

def primes = sieve(Stream from 2, Nil)

然而,这有点难以使用,因为它的每个元素Stream都是使用crosses 列表组成的,列表的数量与素数一样多,并且由于某种原因,似乎这些列表保留在内存中的每个数字Stream.

例如,由注释提示,primes take 6000 contains 56993将抛出GC异常,而primes drop 5000 take 1000 contains 56993在我的测试中会相当快地返回结果.