Prime计数功能的可行实现

Question

任何人都可以提供计算上可行的任何计数功能实现的伪代码吗？我最初尝试编写Hardy-Wright算法,但它的阶乘开始产生可怜的溢出,而其他许多似乎必然会产生类似的问题.我已经在谷歌搜索了实用的解决方案,但最好的是,我发现了非常深奥的数学,这是我在常规程序中从未见过的.

Answer 1

计数函数pi(x)计算不超过x的素数的数量,并且几个世纪以来一直着迷于数学家.在十八世纪初,Adrien-Marie Legendre使用辅助函数phi(x,a)给出了一个公式,该函数计算了不大于x的数字,这些数字没有被第一个素数的筛选所击中; 例如,对于数字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,勒让德的计数公式为pi(x)= phi(x,a)+ a - 1,其中a = pi(floor(sqrt(x))).勒让德使用他的公式来计算pi(10 ^ 6),但是他报告了78526而不是78498的正确答案,尽管错误,但是对于错综复杂的手动计算来说,这是非常接近的.

在20世纪50年代,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 log n和n(log n + log log n)之间的第n个素数P(n)限制为n> 5,因此计算pi在边界并使用二分来确定第n个素数,当边界接近时切换到筛选.

我在博客上讨论了几个条目中的素数.