杜教筛学习笔记

前言

这里假设你已经会了以下东西：

• 狄利克雷卷积
• 莫比乌斯反演
• 数论分块

正文

一切的开始

先来看一个题目（Bzoj 4805）：

给你一个整数$n$，求以下东西：

求$sum_{i=1}^n varphi_i$

可以线性筛吧...但是如果$nleq 2e9$呢？

杜教筛是啥？

设$mathbf S(n)=sum_{i=1}^nmathbf f(i)$，其中$mathbf f$是积性函数。

我们不妨设一个$mathbf g$，不管它具体是啥，做一个狄利克雷卷积：

$$ (mathbf gast mathbf f)(i)=sum_{d|i}mathbf g(d)mathbf f(frac id) $$

然后再求一个和：

$$ sum_{i=1}^nsum_{d|i}mathbf g(d)mathbf f(frac id) \=sum_{d=1}^nmathbf g(d)sum_{d|i}mathbf f(frac id) \=sum_{d=1}^nmathbf g(d)sum_{i=1}^{n/d}mathbf f(i) \=sum_{d=1}^nmathbf g(d)mathbf S(frac nd) $$

然后还有：

$$ mathbf g(1)mathbf S(n)=sum_{i=1}^nmathbf g(i)mathbf S(frac ni)- sum_{i=2}^nmathbf g(i)mathbf S(frac ni) $$

前面减数就是我们之前推的那个东西，后面被减数可以数论分块递归算：

$$ mathbf g(1)mathbf S(n)=sum_{i=1}^n(mathbf gast mathbf f)(i)- sum_{i=2}^nmathbf g(i)mathbf S(frac ni) $$

比如要求$sum_{i=1}^nmu(i)$

$$ mathbf g(1)mathbf S(n)=sum_{i=1}^n(mathbf gast mu)(i)- sum_{i=2}^nmathbf g(i)mathbf S(frac ni) $$

找一个积性函数让他们的狄利克雷卷积很好算：

我们知道$mu$的逆为$1$，也就是说$(1astμ)=epsilon$

那么单位元的前缀和是1，舒服啊，所以取$mathbf g(n)$为$mathbf 1(n)=1$，

所以原式可化为：

$$ mathbf S(n)=1-sum_{i=2}^nmathbf S(frac ni) $$

那么回归到例题

$$ mathbf S(n)=sum_{i=1}^n varphi_i \mathbf g(1)mathbf S(n)=sum_{i=1}^n(mathbf gast varphi)(i)- sum_{i=2}^nmathbf g(i)mathbf S(frac ni) $$

这时我们要找一个$mathbf g$使得$mathbf g$和$mathbf g ast varphi$的前缀和都很好算

这时我们想到了$varphi=mu*mathbf{id}$，所以啊：$varphi*mathbf 1=mathbf{id}$，所以同样取$mathbf g(n)$为$mathbf 1(n)=1$

$$ \mathbf S(n)=frac{n imes (n+1)}{2}- sum_{i=2}^nmathbf S(frac ni) $$

所以，杜教筛的重点就在于选取一个有利于计算的$mathbf g$...

下面给出这题的代码：

#include <cstdio>
#include <cstring>
#include <algorithm>
using std::min; using std::max;
using std::swap; using std::sort;
typedef long long ll;

const int N = 2147483647, M = 2e6 + 10, K = 1.3e3 + 10;
ll phi[M], s1[M], s2[K]; bool vis[K], notprime[M];
int _n, prime[M], tot;

template<typename T>
void read(T &x) {
    int flag = 1; x = 0; char ch = getchar();
    while(ch < '0' || ch > '9') { if(ch == '-') flag = -flag; ch = getchar(); }
    while(ch >= '0' && ch <= '9') x = x * 10 + ch - '0', ch = getchar(); x *= flag;
}

ll S(int n) {
	if(n < M) return s1[n];
	int x = _n / n;
	if(vis[x]) return s2[x];
	vis[x] = true; ll &ans = s2[x];
	ans = 1ll * n * (n + 1) / 2;
	for(int i = 2, j; i <= n; i = j + 1)
		j = n / (n / i), ans -= (j - i + 1) * S(n / i);
	return ans;
}

void get_phi() {
	phi[1] = 1;
	for(int i = 2; i < M; ++i) {
		if(!notprime[i]) prime[++tot] = i, phi[i] = i - 1;
		for(int j = 1; j <= tot && i * prime[j] < M; ++j) {
			notprime[i * prime[j]] = true;
			if(!(i % prime[j])) { phi[i * prime[j]] = phi[i] * prime[j]; break; }
			else phi[i * prime[j]] = phi[i] * (prime[j] - 1);
		}
	}
}

int main () {
	read(_n), get_phi();
	for(int i = 1; i < M; ++i) s1[i] = s1[i - 1] + phi[i];
	printf("%lld
", S(_n));
	return 0; 
} 

一些例题：

• Bzoj3512 DZY Loves Math IV
• LuoguP3768 简单的数学题