javascript 寻找快速素数计数函数

fhity93d  于 2022-11-27  发布在  Java
关注(0)|答案(2)|浏览(173)

我需要计算小于或等于某个N的素数的个数,这就是素数计数函数或PI函数。我有一个,但它太慢了:

function PI(x) {
    var primes = 4;
    for (var i = 3; i <= x; i += 2) {
        if (i % 3 === 0 || i % 5 === 0 || i % 7 === 0) continue;
        var r = ~~Math.sqrt(i), p = true;
        for (var j = 2; j <= r; j++) {
            if (i % j === 0) {
                p = false;
                break;
            }
        }
        if (p)
            primes++;
    }
    return primes;
}

"计算PI(1000000)几乎需要1秒,而计算PI(10000000)对我来说需要20秒,太慢了,有没有更快的方法?"

vfhzx4xs

vfhzx4xs1#

我甚至后来的党,但得到了触发的“素数计数功能”的问题标题,因为这并没有得到太多的关注SO。

TLDR; LMO是如何工作的?跳到结尾查看代码和结果...

素数计数函数技术的主要进步历史如下:
1.大约在1830年,阿德里安-玛丽·勒让德提出了一个公式,使用“包含-排除原理”来计算小于x的素数个数,而不需要列举所有的素数,只需要知道x的平方根的基素数值。

  1. 1870年,E. D. F. Meissel改进了Legendre公式,减少了一些运算次数,但要求知道x^(2/3)以下的素数而不是x^(1/2),这里“知道素数”一般是由SoE来完成的,但只要求存储O(x^(1/3)/(ln x))而不是O(x^(1/2)/(ln x))。虽然它只是将时间复杂度从O(x/((ln x)^2))降低到O(x/((ln x)^3))
    1.最后Lagarias, Miller, and Odzlinko (LMO) developed their algorithm in 1985改进了Meissel方法,仍然需要相同的“素值知识”,直到x^(2/3)的基本版本,但通过将计算分成两组计算,大大减少了计算“Phi”所需的时间(至少,更多的是为了以后的改进),一是决定“普通”树叶的“披”函数,但是这个工作可以忽略,因为只有大约x^(1/3)个这样的函数,并且通过使用部分筛选和计数技术大大减少了计算剩余的“特殊”叶的时间。(x,a)”是不能被任何小于ath素数的素数整除的小于x的所有数的计数,这正是通过筛选基素数值的倍数直到ath素数而产生的结果。基本LMO方法的计算复杂度是O(x^(2/3)(log (log x)))
    为了理解这些技巧是如何使用的,让我们在计算素数到64的平凡范围时,依次考虑它们,如下所示:
    1.考虑素数高达64,我们可以很容易地确定为2,3,5,7,11,13,17,19,23,29,31,37,41,43,47,53,57,和61,总共18个素数;在本练习中,考虑我们只知道对于勒让德技术(2、3、5和7的四个素数),64的平方根为8,或者对于基于迈塞尔的技术,64^(2/3)的范围为16(2、3、5、7、11和13的六个素数)。
    1.勒让德的方法说,我们只需要计算“Phi(64,4)+ 4 - 1”的值,其中4是平方根以下的素数的数目;“Phi”函数是直到该极限64的未排除值的数目,其被包括是因为它们不是小于极限的平方根的四个基素数值的倍数,并且通过如在先前答案中使用的递归技术来计算,在此情况下表示为64 - 64/2 - 64/3 + 64/2/3 - 64/5 + 64/5/2 + 64/5/3 - 64/5/2/3 - 64/7 + 64/7/2 + 64/7/3 - 64/7/2/3 + 64/7/5 - 64/7/5/2 - 64/7/5/3 + 64/7/5/2/3 = 15,当加4减1时给出正确答案18。这被称为“包含-排除原理,”因为具有奇数个素因子的项被排除/减去,并且包含/添加了具有偶数个素因子的项。您可以看到,这需要大量的除法运算,尽管有一些技术,如前一个答案中使用的缓存,可以减少如此小范围内的除法运算,但不幸的是,只能达到一个极限。
    1.对于Meissel技术,计算Phi(64, 2) + 2 - 1(其中2是上至极限4的立方根的素数的数目)要容易得多,因为64 - 64/2 - 64/3 + 64/2/3 = 21在加到2和减到1时是22,但这不是最终答案,因为我们需要减去“P2”项,其是未排除的部分,因此基于素数5和7,它的计算方法是,把64/5 - 3 + 1加上64/7 - 4 + 1的素数的数目加起来,即4。所以,最终的答案是22减去这4对,应该是18。
  2. LMO技术的工作方式与上述完全相同,只是对于这个微不足道的限制,存在一个“特殊”叶,即64/2/3叶,因为2乘以3等于6超过了64的立方根,即4。如果将其实现为LMO,求出商Phi(64/6, 0)的“Phi”,我们不需要因为除数的最小素因子是2,所以根本不进行筛选,并且只计算剩余值的数量,直到商,以获得相同的答案对于叶的“Phi”来说,10。对于这样一个微不足道的限制,“特殊”叶的数量很小,但重要的是,特殊叶的数量仅增加为O((limit^(2/3))/(log (limit^(1/3)))^2)。由于进行筛选的工作是O(limit^(2/3)(log log limit)),这是更大的,这是整个计算复杂度的控制因素。
    那么我们如何实现基本的LMO算法呢?有以下几个步骤:
    1.我们制作了一个表,表中的素数一直到极限的立方根(大小约为(limit^(1/3)) / (ln limit)),这需要的时间可以忽略不计。

1.我们可以使用这些素数通过递归除法技术来计算“Phi”的“S1”普通部分,这花费的时间可以忽略不计,因为只有O(limit^(1/3))中的计算是如此之少。
1.我们从基素数表(最大的表为O(limit^(1/3))大小)中制作一个素数计数值的表,该表的值一直到立方根;我们还制作了一个中间“Phi”计数的表,用于每个部分筛选步骤,直到立方根基素值,以及一个确定的“S2”根叶参数的表,当“S1”值在上面被计算时,其实际上可以被确定;所有这些都花费可忽略的时间,并且最大的(例如最后一个)正好与范围的立方根成比例,这就是为什么在X1 M25 N1 X中使用空间。
1.我们通过页面片段筛选直到limit^(2/3)范围,其中每个页面在高于根“S2”叶表中找到的最小值的时间被一个基素数值部分剔除,累积未剔除值的数量的计数的和直到每个“S2”叶商值。
1.一个基本的页分段的仅奇数位压缩SoE需要大约80秒来计数素数,直到53位限制的2/3次方,但是我们需要使用一个修改的版本,为一些剔除操作的每个新剔除的值更新计数器表,所以它需要大约一半的时间,即大约120秒。
1.然后,在该实现中,对直到每个“S2”点的未剔除值的计数使用批计数器阵列,对于小于由快速“弹出计数”技术计数的批大小的增量的微小最终范围。原始LMO论文描述了使用二叉索引树来实现这一点,但不试图减小树深度,这是它比该实现稍微慢的部分原因,因为由于该× 1 M28 N1 x成本,它将× 1 M27 N1 x因子添加到剔除时间,剔除时更新计数器的因素,这将使剔除时间慢很多倍;对“特殊”叶片的数量的计数和累加对于该实现的上述最大范围仅花费大约几十秒。
1.当对于给定页段已经完成了“S2”计算所需的达到每个单独极限的所有基素数值时,然后将页面段完全剔除到小于该段中表示的最大值和“P2”的平方的基素数值。在该段的范围内的计数值从该完成的筛分页确定,而不必再次筛分相同的范围。
1.当达到range^(2/3)极限的所有页段都完成时,累加器将包含“S2”和“P2”的所需值,使得S1 + S2 - P2 + number of base primes to the cube root of the limit - 1可作为对该limit的素数计数的所需答案而输出。

以下代码实现了基本LMO算法(没有“alpha”优化,这将大大增加代码复杂性),“S2”特殊叶的代码已被松散但不精确地转换为from the "primecount" C++ code by Kim Walisch

第一个
在运行Google Chrome 96版的英特尔Skylake i5-6500 CPU(3.6 GHz,单线程加速率)上,上述代码可以计算范围内的素数,如下所示:
| 限制了|计数器|时间(秒)|
| - -|- -|- -|
| 1e9|小行星50847534|不超过0.006|
| 1e10|小行星455052511| 0.016分|
| 1e11|小行星4118054813| 0.055分|
| 1 e12年|3760791年2018月号|零点二三四|
| 1e13|小行星346065536839|一点一一七|
| 14日|小行星3204941750802|五点五九五|
| 15日|小行星29844570422669|二十八点四三四|
| 2**53-1(约9 e15)|252252704148404号文件|一百四十块二一九|
请注意,当在node.js下运行时,会有大约0.1秒的额外初始化延迟,但在node.js版本16.6.1下运行时会快一些(可能快3%)。

2izufjch

2izufjch2#

时间复杂度为O(n^2/3),空间复杂度为O(n^2/3)。

function eratosthenesWithPi(n) {
  let array = [], upperLimit = Math.sqrt(n), output = [];
  let pi = [0, 0]

  for (let i = 0; i < n; i++) {
    array.push(true);
  }

  for (let i = 2; i <= upperLimit; i++) {
    if (array[i]) {
      for (var j = i * i; j < n; j += i) {
        array[j] = false;
      }
    }
  }

  let cnt = 0

  for (let i = 2; i < n; i++) {
    if (array[i]) {
      output.push(i);
      cnt++
    }

    pi.push(cnt)
  }

  return {primes: new Uint32Array(output), pi: new Uint32Array(pi)}
}

const phiMemo = []
let primes = []

function Phi(m, b) {
  if (b === 0)
    return m
  if (m === 0)
    return 0

  if (m >= 800) {
    return Phi(m, b - 1) - Phi(Math.floor(m / primes[b - 1]), b - 1)
  }

  let t = b * 800 + m

  if (!phiMemo[t]) {
    phiMemo[t] = Phi(m, b - 1) - Phi(Math.floor(m / primes[b - 1]), b - 1)
  }

  return phiMemo[t]
}

const smallValues = [0, 0, 1, 2, 2, 3]
let piValues

function primeCountingFunction(x) {
  if (x < 6)
    return smallValues[x]

  let root2 = Math.floor(Math.sqrt(x))
  let root3 = Math.floor(x ** (1/3))

  let top = Math.floor(x / root3) + 1

  if (root2 + 1 >= primes.length) {
    let res = eratosthenesWithPi(top + 2)

    primes = res.primes
    piValues = res.pi
  }

  let a = piValues[root3 + 1], b = piValues[root2 + 1]

  let sum = 0

  for (let i = a; i < b; ++i) {
    let p = primes[i]

    sum += piValues[Math.floor(x / p)] - piValues[p] + 1
  }

  let phi = Phi(x, a)

  return phi + a - 1 - sum
}

console.log(primeCountingFunction(1e8))

在JSFiddle上试用:https://jsfiddle.net/vo0g274f/1/
对我来说,这需要大约31毫秒。我正在研究一个更节省空间的方法atm,一旦我完成了,我会把它贴在这里。

相关问题