素数计数功能的可行实现答案

【问题标题】：Feasible implementation of a Prime Counting Function [closed]素数计数功能的可行实现
【发布时间】：2013-10-04 22:34:14
【问题描述】：

谁能提供任何prime-counting function 实现的计算上可行的伪代码？我最初尝试对Hardy-Wright algorithm 进行编码，但它的阶乘开始产生可悲的溢出，而且许多其他的似乎必然会产生类似的问题。我在 Google 上搜索了实用的解决方案，但充其量只是发现了非常深奥的数学，我从未见过在传统程序中实现过。

【问题讨论】：

这不是真的floor(x/j) * j == x - (x%j)。然后您链接到的公式变为pi(x) = (-1) + SUM{j=3..n}( ((j-2)!) % j ) (?)。接下来使用模乘法（即5! % 7 == (((((2*3)%7)*4)%7)*5)%7）。
@SeverynKozak 他们说这不是编程问题的原因是因为您的问题不包含代码。
至于提供的不仅仅是伪代码，my answer in another SO thread 提供了 JavaScript 中基本 Legarias-Miller-Odzlinko (LMO) 1985 算法/方法的工作代码 sn-p（可在浏览器中运行） . LMO 构成了所有现代实用的一般素数计数函数的基础，并且已经使用了一种变体 (Gourdon, 2001) 来计算高达 1e28 的素数，尽管即使在现代计算机上使用多线程也需要几个月的时间。大量内核和大量内存。

标签： performance primes pseudocode

【解决方案1】：

素数计数函数 pi(x) 计算不超过 x 的素数，几个世纪以来一直让数学家着迷。 18 世纪初，Adrien-Marie Legendre 给出了一个公式，它使用辅助函数 phi(x,a) 计算不大于 x 且未通过第一个 a 素数筛选的数；例如，对于数字 1、7、11、13、17、19、23、29、31、37、41、43、47 和 49，phi(50,3) = 14。 phi 函数可以计算为 phi (x,a) = phi(x,a-1) - phi(x/P(a),a-1)，其中 phi(x,1) 是不超过 x 和 P(a) 的奇数个数是第 a 个素数（从 P(1)=2 开始计数）。

function phi(x, a)
  if (phi(x, a) is in cache)
    return phi(x, a) from cache
  if (a == 1)
    return (x + 1) // 2
  t := phi(x, a-1) - phi(x // p[a], a-1)
  insert phi(x, a) = t in cache
  return t

数组 p 存储小 a 的第 a 个素数，通过筛分计算得出。缓存很重要；没有它，运行时间将是指数级的。给定 phi，Legendre 的素数计数公式是 pi(x) = phi(x,a) + a - 1，其中 a = pi(floor(sqrt(x)))。 Legendre 使用他的公式计算 pi(10^6)，但他报告了 78526 而不是正确答案 78498，即使错误，对于复杂的手动计算来说也非常接近。

在 1950 年代，Derrick H. Lehmer 给出了一种改进的素数计数算法：

function pi(x)
  if (x < limit) return count(primes(x))
  a := pi(root(x, 4)) # fourth root of x
  b := pi(root(x, 2)) # square root of x
  c := pi(root(x, 3)) # cube root of x
  sum := phi(x,a) + (b+a-2) * (b-a+1) / 2
  for i from a+1 to b
    w := x / p[i]
    lim := pi(sqrt(w))
    sum := sum - pi(w)
    if (i <= c)
      for j from i to lim
        sum := sum - pi(w / p[j]) + j - 1
  return sum

例如，pi(10^12) = 37607912018。即使使用这些算法及其现代变体以及速度非常快的计算机，计算大的 pi 值仍然非常繁琐；在撰写本文时，已知的最大值为 pi(10^24) = 18435599767349200867866。

要使用此算法计算第 n 个素数，素数定理的推论将第 n 个素数 P(n) 限制在 n log n 和 n(log n + log log n) 之间，因为 n > 5 ，因此在边界处计算 pi 并使用二分法确定第 n 个素数，当边界接近时切换到筛分。

我在my blog 的几个条目中讨论素数。

【讨论】：

从第一个代码 sn-p 暗示 phi(x,a)=t，但这里存储了 phi(x, a-1)=t（如果后者为真，则会有一个更快的算法来解决这个问题）。我会自己编辑这个，但我们有编辑的愚蠢限制必须超过 6 个字符。
@StrategyThinker：感谢您的更正。固定。
最近计算了pi(10^25)和pi(10^26)。见page 40 here。
您的帖子对我自己的实现非常有帮助。记忆 pi(x) 也是一个很大的加速。我在计算 a、b、c 时遇到了一些舍入问题，所以我使用了 danaj 的舍入：programmingpraxis.com/2011/07/22/…
行：w := n / p[i]，似乎引用了变量“n”，但我在任何地方都看不到这个变量。我错过了什么吗？

【解决方案2】：

维基百科也可能有所帮助。 prime counting 上的文章包含一些提示。对于初学者，我会推荐 Meissel 在“评估 π(x) 的算法”一节中的算法，这是一种不会生成所有素数的最简单算法。

我还发现 Pomerance 和 Crandall "Prime numbers a computational perspective" 的书很有帮助。这本书对素数计数方法进行了详细且易于理解的描述。但是请记住，对于这里的大多数读者来说，该主题的性质有点太高级了。

【讨论】：