杜教筛

狄利克雷卷积

定义两个函数的狄利克雷卷积如下：

\[ \begin{aligned} f(x) * g(x) &= (f * g)(x)\\ &= \sum_{ab=x}{f(a)g(b)}\\ &= \sum_{d \mid x}{f(d)g(\frac{x}{d})} \end{aligned} \]

其中第一个表达式作为定义式，可以在进行和号变换的时候对理解产生很大帮助。第二个式子则更为常见，可以应用在各种反演之中：

比如常见的莫比乌斯反演和欧拉反演（证明见反演）（好吧没有证明）（\(I(x)=1\)）：

\[ \begin{aligned} (\mu * I)(x) &= \sum_{d \mid x}{\mu(x)I(x)}\\ &= \sum_{d \mid x}{\mu(x)}\\ &= [x=1]\\ \\(\varphi * I)(x) &= \sum_{d \mid x}{\varphi(x)I(x)}\\ &=\sum_{d \mid x}{\varphi(x)}\\ &=x \end{aligned} \]

好，这些知识就足够了。

杜教筛

杜教筛其实就是能筛这样的一个函数：

\[ f(x)=G(x)-\sum_{d=2}^{x}f(\lfloor\frac{x}{d}\rfloor) \]

如果 \(G(x)\) 可以很快地算出来的话，那么第二项可以利用数论分块递归 \(\sqrt x\) 次。然后对于这个东西记忆化。时间复杂度我不会证，这里贴一个搜出来的 。

复杂度证明：
我们使用同一法来简单证明。总时间复杂度可以看成是所有在递归过程中出现过的 \(n\) 对应的时间复杂度的和。
显然并不是每一个 \([1, n]\) 之间的数都出现了，具体哪些出现了呢？大致是所有 \(\lfloor\frac{n}{i}\rfloor\) 在每个 \(\lfloor\frac{n}{i}\rfloor\) 上花费的时间又是 \(O(\sqrt{\lfloor\frac{n}{i}\rfloor})\) 的，因此： \[T(n)=\sum_{i=1}^{\sqrt n}{\sqrt{\frac{n}{i}} + \sqrt i}\]
显然这个式子中 \(\sqrt i\) 的项每一项都不大于 \(\sqrt{\frac{n}{i}}\) 的项，因此可以忽略。 那么 \[T(n)=\sum_{i=1}^{n}{\sqrt{\frac{n}{i}}} \le \int_1^{\sqrt n}{\frac{\sqrt n}{\sqrt x}} = 2 \sqrt n \sqrt[4] x = O(n^{\frac{3}{4}})\]

然后我们如果用线性筛求出前 \(O(n^{\frac{2}{3}})\) 个前缀和，就可以进一步将复杂度压到 \(n^{\frac{2}{3}}\) 。

那么问题来了，什么样的函数可以满足 \(f(x)=G(x)-\sum_{d=2}^{x}f(\lfloor\frac{x}{d}\rfloor)\) 呢？

除了题目直接告诉你的性质之外，一般的函数需要借用之前提到的卷积，这样可以求出它的前缀和（记为 \(S(x)\) ）。我们先找到一个函数 \(g(x)\) ，满足 \(g(x)\) 和 \((f \* g)(x)\) 都可以快速算出。比如对于 \(\mu\) 函数和 \(\varphi\) 函数来说， \(I\) 函数就是一个合理的选择。因为 \(I(x)\) 就是 1 ，而 \((\mu \* I)(x)=[x=1],(\varphi*I)(x)=x\) 。

具体怎么实现呢？我们来看：

\[ \begin{aligned} \sum_{i=1}^{n}{(f*g)(i)}&=\sum_{i=1}^{n}{\sum_{d \mid i}{f(\frac{i}{d})g(d)}}\\ \text{变换枚举顺序，}&\text{新的 }i\text{ 枚举原来的 }d\text{ ，新的 }j\text{ 枚举原来的 }\frac{i}{d}\\ \text{可以发现所有的}&\text{二元组 }(i,j)\text{ 就是所有满足 }i*j\le n\text{ 的有序二元组，所以有：}\\ \sum_{i=1}^{n}{(f*g)(i)}&=\sum_{i=1}^{n}{\sum_{j=1}^{\lfloor\frac{n}{i}\rfloor}{g(i)f(j)}}\\ &=\sum_{i=1}^{n}{g(i)\sum_{j=1}^{\lfloor\frac{n}{i}\rfloor}{f(j)}}\\ &=\sum_{i=1}^{n}{g(i)S(\lfloor\frac{n}{i}\rfloor)} \end{aligned} \]

移项变形：

\[ S(n)=\frac{\sum_{i=1}^{n}{(f*g)(i)} - \sum_{i=2}^{n}{g(i)S(\lfloor\frac{n}{i}\rfloor)}}{g(1)} \]

下面给出 \(f=\mu,g=I\) 以及 \(f=\varphi,g=I\) 的特例。推导过程可以给读者作为练习。

\[ \begin{aligned} S_{\mu}(n)&=1 - \sum_{i=2}^{n}{S(\lfloor\frac{n}{i}\rfloor)}\\ S_{\varphi}(n)&=\frac{n(n+1)}{2} - \sum_{i=2}^{n}{S(\lfloor\frac{n}{i}\rfloor)} \end{aligned} \]

好的。我们再来看一些更有意思的函数：求 \(f(x)= x \mu(x)\) 的前缀和。

我们可以令 \(S(x)=\sum_{i=1}^{x}{f(x)},g(x)=x\)

那么有：

\[ \begin{aligned} S(n)&=\frac{\sum_{i=1}^{n}{(f*g)(i)} - \sum_{i=2}^{n}{g(i)S(\lfloor\frac{n}{i}\rfloor)}}{g(1)}\\ &=\frac{\sum_{i=1}^{n}{\sum_{d \mid i}{g(d)f(\frac{i}{d})}} - \sum_{i=2}^{n}{g(i)S(\lfloor\frac{n}{i}\rfloor)}}{g(1)}\\ &=\sum_{i=1}^{n}{\sum_{d \mid i}{d\frac{i}{d}\mu(\frac{i}{d})}} - \sum_{i=2}^{n}{iS(\lfloor\frac{n}{i}\rfloor)}\\ &=\sum_{i=1}^{n}{\sum_{d \mid i}{i\mu(\frac{i}{d})}} - \sum_{i=2}^{n}{iS(\lfloor\frac{n}{i}\rfloor)}\\ &=\sum_{i=1}^{n}{i\sum_{d \mid i}{\mu(d)}} - \sum_{i=2}^{n}{iS(\lfloor\frac{n}{i}\rfloor)}\\ &=\sum_{i=1}^{n}{i[i=1]} - \sum_{i=2}^{n}{iS(\lfloor\frac{n}{i}\rfloor)}\\ &=1 - \sum_{i=2}^{n}{iS(\lfloor\frac{n}{i}\rfloor)} \end{aligned} \]

用同样的方法可以求出 \(\varphi(i)i\) 的前缀和。

例题

P4213【模板】杜教筛（Sum）

AC 代码

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
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71

#include <bits/stdc++.h>
using namespace std;

map<long long, long long> sum_mu, sum_phi;
map<long long, long long>::iterator it;

bool not_prime[1664500];
long long phi[1664500], mu[1664500];
vector<int> primes;

long long calc_sum_mu(long long x) {
	if (x < 1664500) return mu[x];
	it = sum_mu.find(x);
	if (it != sum_mu.end()) return it -> second;
	long long res = 1;
	for (long long i = 2, j; i <= x; i = j + 1) {
		j = x / (x / i);
		res -= (j - i + 1) * calc_sum_mu(x / i);
	}
	return sum_mu[x] = res;
}

long long calc_sum_phi(long long x) {
	if (x < 1664500) return phi[x];
	it = sum_phi.find(x);
	if (it != sum_phi.end()) return it -> second;
	long long res = (x * (x + 1)) >> 1;
	for (long long i = 2, j; i <= x; i = j + 1) {
		j = x / (x / i);
		res -= (j - i + 1) * calc_sum_phi(x / i);
	}
	return sum_phi[x] = res;
}

int main() {
	// freopen(".in", "r", stdin);
	// freopen(".out", "w", stdout);
	mu[1] = 1;
	phi[1] = 1;
	for (int i = 2; i < 1664500; i++) {
	    if (!not_prime[i]) {
			primes.emplace_back(i);
			phi[i] = i - 1;
			mu[i] = -1;
		}
		for (int p : primes) {
			if (i * p >= 1664500) break;
			not_prime[i * p] = true;
			if (i % p == 0) {
				phi[i * p] = phi[i] * p;
				mu[i * p] = 0;
				break;
			}
			phi[i * p] = phi[i] * phi[p];
			mu[i * p] = -mu[i];
		}
	}
	for (int i = 2; i < 1664500; i++) {
		phi[i] += phi[i - 1];
		mu[i] += mu[i - 1];
	}
	int _;
	scanf("%d", &_);
	while (_--) {
		static long long x;
		scanf("%lld", &x);
		printf("%lld %lld\n", calc_sum_phi(x), calc_sum_mu(x));
		// cerr << mu[int(x)] << endl;
	}
	return 0;
}