转载自https://oi-wiki.org/math/du/
在莫比乌斯反演的题目中,往往要求出一些数论函数的前缀和,利用 杜教筛 可以快速求出这些前缀和。
杜教筛
求 $displaystyle S(n)=sum_{i=1}^n f(i)$
我们要想办法构造一个 $S(n)$ 关于 $S(left lfloor frac{n}{i} ight floor)$.
构造两个积性函数 $g$,使得 $h=f*g$ 的前缀和很好求
$$egin{aligned}
sum _{i=1}^n h(i) &= sum_{i=1}^n sum _{d|n} g(d)f(frac{n}{d}) \
&= sum_{d=1}^n g(d)cdot sum_{i=1}^{frac{n}{d}} f(i) \
&= sum_{d=1}^n g(d)S(left lfloor frac{n}{d}
ight
floor) \
&= sum _{d=2}^n g(d)S(left lfloor frac{n}{d}
ight
floor) + g(1)S(n)
end{aligned}$$
因此 $g(1)S(n) = sum_{i=1}^nh(i) - sum _{d=2}^n g(d)S(left lfloor frac{n}{d} ight floor)$.
那么假设我们可以快速对 $sum_{i=1}^n(f*g)(i)$ 求和,并用数论分块求解 $sum_{i=2}^n g(i)S(frac{n}{i})$ 就可以在较短时间内求得 $g(1)S(n)$.
莫比乌斯函数前缀和
前面证明过如下卷积式:$mu * 1 = varepsilon$
$$egin{aligned}
herefore S_1(n) &= sum_{i=1}^n varepsilon (i) - sum_{i=2}^nS_1(left lfloor frac{n}{i}
ight
floor) \
&= 1 - sum_{i=2}^nS_1(left lfloor frac{n}{i}
ight
floor)
end{aligned}$$
观察到 $left lfloor frac{n}{i} ight floor$ 最多只有 $O(sqrt n)$ 种取值,我们就可以应用整除分块(或称数论分块)来计算每一项的值了。
直接计算的时间复杂度为 $O(n^{frac{3}{4}})$。考虑先线性筛预处理出前 $n^{frac{2}{3}}$ 项,剩余部分的时间复杂度为 $O(int _0^{n^frac{1}{3}} sqrt{frac{n}{x}}dx) = O(n^{frac{2}{3}})$.
对于较大值,需要用 map (或者unordered_map)存在其对应的值,方便以后使用时直接使用之前计算的结果。
欧拉函数前缀和
当然也可以用杜教筛求出 $varphi (x)$ 的前缀和,但是更好的方法是应用莫比乌斯反演:
$displaystyle egin{array}
{l}{sum_{i=1}^{n} sum_{j=1}^{n} 1[operatorname{gcd}(i, j)=1]=sum_{i=1}^{n} sum_{j=1}^{n} sum_{d|i, d| j} mu(d)} \
{=sum_{d=1}^{n} mu(d)leftlfloorfrac{n}{d}
ight
floor^{2}}
end{array}$
由于题目所求的是 $sum_{i=1}^{n} sum_{j=1}^{i} 1[operatorname{gcd}(i, j)=1]$,所以我们排除掉 $i=j=1$ 的情况,并将结果除2即可。
观察到,只需求出莫比乌斯函数的前缀和,就可以快速计算出欧拉函数的前缀和了。时间复杂度为 $O(n^{frac{2}{3}})$
使用杜教筛求解
$ecause varphi * 1=I D$
$ herefore S(n)=frac{1}{2} n(n+1)-sum_{i=2}^{n} Sleft(leftlfloorfrac{n}{i} ight floor ight)$
代码实现
//求莫比乌斯函数和欧拉函数的前缀和
#include <algorithm> #include <cstdio> #include <cstring> #include <map> using namespace std; const int maxn = 2000010; typedef long long ll; ll T, n, pri[maxn], tot, mu[maxn], sum_mu[maxn]; bool vis[maxn]; map<ll, ll> mp_mu; //可换成unordered_map ll S_mu(ll x) { if (x < maxn) return sum_mu[x]; if (mp_mu[x]) return mp_mu[x]; ll ret = 1ll; for (ll i = 2, j; i <= x; i = j + 1) { j = x / (x / i); ret -= S_mu(x / i) * (j - i + 1); } return mp_mu[x] = ret; } ll S_phi(ll x) { ll ret = 0ll; for (ll i = 1, j; i <= x; i = j + 1) { j = x / (x / i); ret += (S_mu(j) - S_mu(i - 1)) * (x / i) * (x / i); } return ((ret - 1) >> 1) + 1; } void initMu() { mu[1] = 1; for (int i = 2; i < maxn; i++) { if (!vis[i]) pri[++tot] = i, mu[i] = -1; for (int j = 1; j <= tot && i * pri[j] < maxn; j++) { vis[i * pri[j]] = true; if (i % pri[j] == 0) { mu[i * pri[j]] = 0; break; } else { mu[i * pri[j]] = -mu[i]; } } } for (int i = 1; i < maxn; i++) sum_mu[i] = sum_mu[i - 1] + mu[i]; } int main() { initMu(); scanf("%lld", &T); while (T--) { scanf("%lld", &n); printf("%lld %lld ", S_phi(n), S_mu(n)); } return 0; }