题意

输出[1..n]的质数个数 (1 <= n <= 1e11) 。时间限制6s,空间限制64M 。

题解

很显然,要用线性筛的话时间和空间都不够。

有一种Meissel–Lehmer算法可以解决该问题,不过本文要介绍的是一种dp解法。

定义$SR(n,p)$为,$2..n$被小于等于$p$的质数筛后剩下的数的个数;也就是说,在n的范围内,是质数,或者是只由大于$p$的质数相乘得到的数的个数。
记小于等于$n$的质数个数为$\pi(n)$,那么显然有$SR(n,n)=\pi(n)$ 。

对$SR(n,p)$分两类讨论:

  • 当$p$不是质数或者$n<p^2$时,有$SR(n,p)=SR(n,p-1)$;
  • 当$p$是质数且$n\ge p^2$时,$SR(n,p)$ 也可由 $SR(n,p-1)$ 推得:$$SR(n,p)= SR(n,p−1)− SR(\frac{n}{p}, p−1)+ SR(p−1,p−1)$$ 表示要从 $2..p-1$ 筛后剩下的数中去掉那些质因子均大于等于 $p$ 且含 $p$ 的数,因为若有小于 $p$ 的质因子则该数已经被筛去了,同理,该数除以p也一定不会被 $2..p-1$ 筛去,所以减去 $SR(\frac{n}{p}, p−1)$,与此同时 $2..p-1$ 的质数也被减掉了,所以加上 $SR(p-1,p-1)$ 。

当然因为空间限制肯定不能直接这样dp,考虑到整个过程中只用到了 $p$ 和 $\frac{n}{p}$,我们直接分为两段dp, 用H[k]表示 $k\le \sqrt{n}$ 时 $SR(\frac{n}{k},p)$ 的值,用L[k]表示 $k \le \sqrt{n}$ 时 $SR(k,p)$ 的值。

先放出代码,原作者*adkroxx*:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
#include <bits/stdc++.h>
using namespace std;

typedef long long ll;
const int N=320000;
ll H[N],L[N];

ll SR(ll n,ll p)
{
ll m;
for (m=1;m*m<=n;++m) {
H[m]=n/m-1;
}
for (ll i=1;i<=m;++i) {
L[i]=i-1;
}
for (ll i=2;i<=m;++i){
if (L[i]==L[i-1]) continue;
for (ll j=1;j<=min(m-1,n/i/i);++j) {
if(i*j<m) H[j]-=H[i*j]-L[i-1];
else H[j]-=L[n/i/j]-L[i-1];
}
for (ll j=m;j>=i*i;--j) {
L[j]-=L[j/i]-L[i-1];
}
}
return H[1];
}

int main()
{
ll n;
while (scanf("%lld",&n)!=EOF) {
printf("%lld\n",SR(n,n));
}
}

L[],H[]初始化为筛内总的个数(不含1)。代码中的i用来枚举质数(L[i]==L[i-1]表示i不是质数),然后对于 $j>i^2$ 的 $j$ 都更新一遍dp值,只是放在了两个不同的数组而已:对于 $j\le \sqrt{n}$ 的$j$更新L[],对于 $j>\sqrt{n}$ 的 $j$ 更新H[]。只要理解了递推式代码还是很好理解的。

这个算法的时间复杂度为$O(n^{\frac{3}{4}})$(证明可以看参考资料内的下一条评论),虽然不如Meissel–Lehmer算法,但是dp的思想非常巧妙,值得学习。

后记

看到csdn和cnblogs上很多人都是把这份代码当作模板直接贴出来不讲解,不仅不尊重原作者,也有悖大家学习数据结构与算法的初衷,是以我不仅要翻译,还要加上一些自己的理解~~

参考资料

Editorial of Educational Codeforces Round 12