杜教筛

杜教筛可以在 $O(n^\frac23)$ 的时间复杂度内求出部分数论函数的前缀和。

假设我们要求出 $S(n) = \sum_{i = 1}^n f(i)$，我们需要尝试构造出一个合适的数论函数 $g$，那么有：

$\begin{aligned} \sum_{i = 1}^n (f * g)(i) &= \sum_{i = 1}^n \sum_{j | i} f\left(\frac ij\right)g(j)\\ &= \sum_{j = 1}^n g(j) \sum_{i = 1}^{\lfloor\frac nj\rfloor} f(i)\\ &= \sum_{i = 1}^n g(i)S\left(\left\lfloor \frac ni \right\rfloor\right) \end{aligned}$

那么有：

$g(1)S(n) = \sum_{i = 1}^n (f * g)(i) - \sum_{i = 2}^n g(i)S\left(\left\lfloor \frac ni \right\rfloor\right)\\ S(n) = \frac{\sum_{i = 1}^n (f * g)(i) - \sum_{i = 2}^n g(i)S\left(\left\lfloor \frac ni \right\rfloor\right)}{g(1)}\\$

如果 $g, f * g$ 都能快速求前缀和，我们就可以使用数论分块快速递推出 $S(n)$。

直接递推时间复杂度是 $O(n^\frac 34)$ 的，但是如果我们预处理出前 $n^\frac 23$ 的 $S(i)$，那么这个东西的时间复杂度会变成 $O(n^\frac 23)$。

运用

【模板】杜教筛

这题中 $f$ 分别为 $\mu, \varphi$，我们一个一个来。

如果 $f = \mu$，那么我们想到 $\mu * \mathrm{Id}_0 = \epsilon$，而后两者的前缀和都可以 $O(1)$ 求，那么我们令 $g = \mathrm{Id}_0$ 即可。

如果 $f = \varphi$，那么我们想到 $\varphi * \mathrm{Id}_0 = \mathrm{Id}_1$，而后两者的前缀和都可以 $O(1)$ 求，那么我们也令 $g = \mathrm{Id}_0$ 即可。

然后就做完了。

#include <algorithm>
#include <cstring>
#include <cstdio>
#include <vector>
#include <cmath>
using namespace std;
template<typename T, const int M>
struct BHasher{inline int operator()(T x){return (x % M + M) % M;}};
template<typename T1, typename T2, int M, typename Hsr = BHasher<T1, M>>
struct HashMap{
	struct node{int nx;T1 ky;T2 vl;};
	int h[M];vector<node> nd;Hsr hsh;
	HashMap(){memset(h, -1, sizeof(h));}
	void clear(){
		for (auto oo : nd) h[hsh(oo.ky)] = -1;
		nd.clear();
	}
	inline bool find(T1 k){int p = hsh(k);
		for (int i = h[p];i != -1;i = nd[i].nx)
		if (nd[i].ky == k) return 1;return 0;
	}
	inline T2& operator[](T1 k){int p = hsh(k);
		for (int i = h[p];i != -1;i = nd[i].nx)
		if (nd[i].ky == k) return nd[i].vl;
		nd.push_back((node){h[p], k, T2()});
		return nd[h[p] = nd.size() - 1].vl;
	}
};
typedef long long ll;
const int N = 2000005;bool isnp[N];ll sph[N];
int t, n, m, cp, pr[N], mu[N], phi[N], smu[N];
HashMap<int, ll, 1145143> Sp;
HashMap<int, int, 1145143> Sm;
ll Sph(int U){
	if (U <= m) return sph[U];
	if (Sp.find(U)) return Sp[U];
	ll ans = U * (U + 1ll) >> 1;
	for (int l = 2, r;l > 0 && l <= U;l = r + 1){
		r = U / (U / l), ans -= (r - l + 1) * Sph(U / l);
		if (r == U) break;
	}return Sp[U] = ans;
}
int Smu(int U){
	if (U <= m) return smu[U];
	if (Sm.find(U)) return Sm[U];
	int ans = 1;
	for (int l = 2, r;l > 0 && l <= U;l = r + 1){
		r = U / (U / l), ans -= (r - l + 1) * Smu(U / l);
		if (r == U) break;
	}return Sm[U] = ans;
}
int main(){scanf("%d", &t);
	while (t --){scanf("%d", &n), m = floor(pow(n, 0.6666));
		memset(isnp, 0, sizeof(isnp)), isnp[0] = isnp[1] = phi[1] = mu[1] = 1, cp = 0;
		for (int i = 2;i <= m;i ++){
			if (!isnp[i]) pr[++ cp] = i, mu[i] = -1, phi[i] = i - 1;
			for (int j = 1;j <= cp && i * pr[j] <= m;j ++){
				isnp[i * pr[j]] = 1;
				if (i % pr[j] == 0){mu[i * pr[j]] = 0, phi[i * pr[j]] = phi[i] * pr[j];break;}
				mu[i * pr[j]] = -mu[i], phi[i * pr[j]] = phi[i] * phi[pr[j]];
			}
		}
		for (int i = 1;i <= m;i ++) smu[i] = smu[i - 1] + mu[i], sph[i] = sph[i - 1] + phi[i];
		printf("%lld %d\n", Sph(n), Smu(n));
	}
}