从零开始学杜教筛

实际上略过了一小部分基础数论。
杜教筛时间复杂度的证明还需要一点微积分基础。

目录

• 积性函数
  • 定义
• 狄利克雷卷积
  • 定义
  • 性质
  • 狄利克雷逆
  • 应用
  • 恒等式
• 数论分块（整除分块）
  • 重要观察
  • 枚举取值
• 线性筛
  • 原理
  • 正确性
  • 复杂度
  • 代码实现
• 杜教筛
  • 简介
  • 式子推导
  • 时间复杂度
  • 应用
    • 筛 \(\mu\)
    • 筛 \(\varphi\)
  • 代码实现

积性函数

数论函数是定义域在整数上的函数，其中积性函数因为其具有十分朴实且优美的性质，成为了 OI 数论中不可或缺的基础。

定义

设数论函数 \(f(n)\)，若对于任意 \(a\perp b\) 满足 \(f(a)\cdot f(b)=f(ab)\) 则称 \(f(n)\) 为积性函数。
特别的，若对于任意 \(a,b\) 满足 \(f(a)\cdot f(b)=f(ab)\)，则称 \(f(n)\) 为完全积性函数。

常见的积性函数有：

• 单位函数：\(\varepsilon(n)=[n=1]\)（完全积性）
• 恒等函数：\(\text{id}_k(n)=n^k\)（完全积性），当 \(k=1\) 简记为 \(\text{id}(n)=n\)
• 常数函数：\(1(n)=1\)（完全积性）
• 除数函数：\(\sigma_k(n)=\sum_{d|n}d^k\)，\(\sigma_0(n)\) 即因数个数，简记为 \(\text d(n)\)；\(\sigma_1(n)\) 即因子之和，简记为 \(\sigma(n)\)
• 欧拉函数：\(\varphi(n)=\sum_{d=1}^n[gcd(d,n)=1]\)，即 \(n\) 以内与 \(n\) 互质的数的个数。
• 莫比乌斯函数： \(\mu(x)=\begin{cases} 1&(x=1) \\ (-1)^k&(x没有平方数因子，且x的质因子个数为k)\\ 0&(x有平方数因子) \end{cases}\)

同时有两个重要的推论：

1. 积性函数 \(f\) 一定满足 \(f(1)=1\)，因为 \(f(1)=f(1\times1)=f(1)\times f(1)\)。
2. 通过所有质数处的点值可以唯一确定完全积性函数；通过全部 \(p^k\) 处点值可以唯一确定积性函数，由于唯一分解定理。

狄利克雷卷积

定义

对于两个数论函数 \(f,g\)，他们的狄利克雷卷积定义为

\[(f\otimes g(n))=\sum_{d|n}f(d)g({n\over d}) \]

性质

狄利克雷卷积有很多优良的性质：

• 交换律：\(f\otimes g= g\otimes f\)
• 结合律：\(f\otimes g\otimes h=f\otimes(g\otimes h)\)
• 单位元：取 \(\varepsilon(n)=[n=1]\)，则对任意数论函数
  \(f\) 都有 \(f\otimes \varepsilon=\varepsilon\otimes f=f\)
• 两个积性函数的狄利克雷卷积仍然是积性函数。

狄利克雷逆

已知数论函数 \(f\)，求一个数论函数 \(g\) 使 \(f\otimes g=\varepsilon\)。

应用

这个东西有啥用呢？我们不妨试着找找 \(1(n)=1\) 的狄利克雷逆。
记 \(\mu=1^{-1}\)，我们简单计算一下前几项：

\(n\)	\(1\)	\(2\)	\(3\)	\(4\)	\(5\)	\(6\)	\(7\)	\(8\)	\(9\)	\(10\)	\(11\)	\(12\)	\(13\)	\(14\)	\(15\)	\(16\)
\(\mu(n)\)	\(1\)	\(-1\)	\(-1\)	\(0\)	\(-1\)	\(1\)	\(-1\)	\(0\)	\(0\)	\(\ 1\)	\(-1\)	\(\ 0\)	\(-1\)	\(\ 1\)	\(\ 1\)	\(\ 0\)

设 \(p,q\) 为互异质数，我们可以简单观察到一些规律：

• 在 \(n=p\) 处 \(\mu(p)=-1\)
• 在 \(n=pq\) 处 \(\mu(pq)=1\)
• 在 \(n=p^k\ (k>1)\) 处 \(\mu(p^k)=0\)

这些条件看上去都是必要的，但是并没有充分性。
根据狄利克雷卷积的性质，由于 \(1\) 和 \(\varepsilon\) 都是积性函数，所以 \(\mu\) 也是。并且由积性函数推论 2，这些 \(p^k\) 处的点值已经唯一确定了一个积性函数。归纳得出：

\[\mu(x)=\begin{cases} 1&(x=1) \\ (-1)^k&(x没有平方数因子，且x的质因子个数为k)\\ 0&(x有平方数因子) \end{cases}\]

事实上，\(1(n)\) 的狄利克雷逆就是 \(\mu(n)\) 的定义。

我们考虑另一个狄利克雷卷积：\(\mu\otimes\text{id}\)，其中 \(\text{id}(n)=n\)。
记 \(\varphi=\mu\otimes\text{id}\)，我们简单计算一下前几项：

\(n\)	\(1\)	\(2\)	\(3\)	\(4\)	\(5\)	\(6\)	\(7\)	\(8\)	\(9\)	\(10\)	\(11\)	\(12\)	\(13\)	\(14\)	\(15\)	\(16\)
\(\varphi(n)\)	\(1\)	\(1\)	\(2\)	\(2\)	\(4\)	\(2\)	\(6\)	\(4\)	\(6\)	\(\ 4\)	\(10\)	\(\ 4\)	\(12\)	\(\ 6\)	\(\ 8\)	\(\ 8\)

设 \(p\) 为质数，我们可以简单观察到一些规律：

• 在 \(n=p\) 处 \(\mu(p)=p-1\)
• 在 \(n=p^k\ (k>1)\) 处 \(\mu(p^k)=(p-1)p^{k-1}\)

归纳得出：

\[\varphi(n)=n\prod_{i=1}^k{p_i-1\over p_i},\ \text {where}\ n=\prod_{i=1}^kp_i^{r_i} \]

事实上，\((\mu\otimes\text{id})(n)\) 就是 \(\varphi(n)\) 的定义。

恒等式

我们有以下狄利克雷卷积恒等式：

• \(\mu\otimes 1=\varepsilon\)
• \(\mu\otimes\text{id}=\varphi\)
• \(\varphi\otimes1=\text{id}\)
• \(1\otimes1=d=\sigma_0\)
• \(1\otimes\text{id}_k=\sigma_k\)

其中前三个恒等式非常重要，在杜教筛式子的推导中也会出现。

数论分块（整除分块）

这是几乎所有数论题都会用到的计算方法，非常重要。一切来源于一个精巧的观察。

重要观察

对一个正整数 \(n\)，所有 \(\lfloor{n\over i}\rfloor\) 的取值有 \(O(\sqrt n)\) 种。

证明：当 \(i\le\sqrt x\) 时，一一对应 \(O(\sqrt n)\) 种不同的取值；当 \(i>\sqrt n\) 时，由于 \(\lfloor{n\over i}\rfloor\le\sqrt n\)，所以此时 \(\lfloor{n\over i}\rfloor\) 也只有 \(O(\sqrt n)\) 种取值。

枚举取值

现在我们有一个 \(l\) 表示一段取值相同的数的左端点，要求这段数的右端点 \(r\)。即已知 \(\lfloor{n\over l}\rfloor=\lfloor{n\over r}\rfloor\)，求 \(r_{max}\)。

设 \(k=\lfloor{n\over l}\rfloor=\lfloor{n\over r}\rfloor\)，有 \(k\le{n\over r}<k+1\)。同时取倒数，得到 \({1\over k+1}<{r\over n}\le {1\over k}\)。求 \(r\) 上界拿出右边不等式，化成 \(r\le{n\over k}\)。由于求的 \(r\) 是一个整数，不妨加个下取整得到 \(r\le\lfloor{n\over k}\rfloor=\lfloor{n\over \lfloor{n\over l}\rfloor}\rfloor\)。即 \(r_{max}=\lfloor{n\over \lfloor{n\over l}\rfloor}\rfloor\)。

我们可以不断地求出 \(r=\lfloor{n\over \lfloor{n\over l}\rfloor}\rfloor\)，统计整段贡献，再令 \(l = r + 1\)。这样我们就做到了 \(O(\sqrt n)\) 快速求值。

线性筛

原理

顾名思义，线性筛可以在严格线性时间复杂度筛质数或各种积性函数。其原理是用每个数的最小质因子来筛去它。

正确性

具体实现是维护从小到大的质数 \(p_j\)，每次筛去当前数 \(i\times p_j\)，此时默认了 \(p_j\) 是 \(i\times p_j\) 的最小质因子。当枚举到 \(p_j\) 整除 \(i\) 时，我们发现：如果之前没有枚举到 \(i\) 的其他质因子，那么 \(p_j\) 是 \(i\) 的最小质因子；在这之后，如果有另一个更大的 \(p_{j'}\) 整除 \(i\)，那么 \(i\times p_{j'}\) 的最小质因子仍然是 \(p_j\) 而不是 \(p_{j'}\)，我们不能用 \(p_{j'}\) 来筛 \(i\times p_{j'}\)，这样不符合原理。所以在第一处 \(p_j\) 整除 \(i\) 时就需要终止对于 \(i\times p_j\) 的筛。而在此之前筛去的 \(i\times p_{j''}\) 如果有更小的质因子，那一定是 \(i\) 更小的质因子，在之前就会终止。所以我们的默认也是正确的。

复杂度

复杂度如何保证？对于每个数去筛都有一次终止，即 \(O(n)\) 次终止；在终止之前，每一个数 \(i\times p_j\) 都只会被最小质因子筛一次，同时这个 \(i\times p_j\) 只会被最小质因子打上一次非质数的标记。这几个行为都是 \(O(n)\) 的，故总复杂度也是 \(O(n)\)。

代码实现

贴一份线性筛代码，其实很简短。

线性筛 $\mu$ 和 $\varphi$

void init() {
	mu[1] = 1, phi[1] = 1;
	for(int i = 2; i <= maxv; i++) {
		if(isp[i]) {
			primes.push_back(i);
			mu[i] = -1, phi[i] = i - 1;
		}
		for(int p : primes) {
			if(i * p > maxv) break;
			isp[i * p] = false;
			if(i % p) {
				mu[i * p] = -mu[i], phi[i * p] = phi[i] * (p - 1);
			}
			else {
				mu[i * p] = 0, phi[i * p] = phi[i] * p;
				break;
			}
		}
	}
	return;
}

杜教筛

简介

对于数论函数 \(f\)，杜教筛可以在低于线性时间复杂度（\(O(n^{3\over4})\) 甚至 \(O(n^{2\over3})\)）计算 \(S(n)=\sum_{i=1}^nf(i)\) 的点值。

式子推导

找另一个恰当的数论函数 \(g\)（通常是积性函数），考虑 \(f\otimes g\)：

\[\sum_{i=1}^n(f\otimes g)(i)=\sum_{i=1}^n\sum_{d|i}f(d)g({i\over d}) \]

为了避免枚举因数，我们更换指标：用新的 \(d'\) 代替 \(i\over d\)，新的 \(i'\) 代替 \(d\)。由于 \(d'\) 是 \(i\) 的因数，如果枚举 \(d'\)，所有合法的 \(i'\) 表示的是 \(d'\) 在 \(n\) 以内合法的倍数，共有 \(\lfloor{n\over d'}\rfloor\) 个。于是：

\[\sum_{i=1}^n\sum_{d|i}f(d)g({i\over d})=\sum_{d'=1}^ng(d')\sum_{i'=1}^{\lfloor{n\over d'}\rfloor}f(i') \]

把对 \(f(i')\) 求和换成 \(S(\lfloor{n\over d'}\rfloor)\)，再把 \(d'\) 换成更简单的指标：

\[\sum_{i=1}^ng(i)S(\lfloor {n\over i}\rfloor) \]

（以上部分也可以用于在其他时候处理狄利克雷卷积）

接下来我们处理 \(S(n)\)。如果你选取的 \(g\) 是积性函数，那么 \(S(n)=g(1)S(n)\) 可以写成 \(\sum_{i=1}^ng(i)S(\lfloor{n\over i}\rfloor)-\sum_{i=2}^ng(i)S(\lfloor{n\over i}\rfloor)\)（如果不是，则整体除以 \(g(1)\)）。即有：

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

令 \(f\otimes g=h\)，则：

\[g(1)S(n)=\sum_{i=1}^nh(i)-\sum_{i=2}^ng(i)S(\lfloor{n\over i}\rfloor) \]

如果 \(h\) 和 \(g\) 的前缀和好求，那么我们只需要数论分块就可以快速求出 \(S(\lfloor{n\over i}\rfloor)\)。所以一个恰当的 \(g\) 条件是 \(f\otimes g\) 和 \(g\) 的前缀和都好求。

时间复杂度

令 \(R(n)=\{\lfloor{n\over k}\rfloor:k=2,3,\cdots,n\}\)

我们认为求 \(h\) 和 \(g\) 的前缀和是 \(O(1)\) 的。每个 \(S(k)\ (k\in R(n))\) 均只会计算一次，其中 \(|R(n)|=O(\sqrt n)\)（数论分块结论）。

设计算一次 \(S(n)\) 的时间复杂度为 \(T(n)\)，则：

\[\begin{aligned} T(n)&=\sum_{k\in R(n)}T(k)\\ &=\Theta(\sqrt n)+\sum_{k=1}^{\lfloor\sqrt n\rfloor}O(\sqrt k)+\sum_{k=2}^{\lfloor\sqrt n\rfloor}O(\sqrt {n\over k})\\ &=O(\int_0^{\sqrt n}(\sqrt x+\sqrt{n\over x})\text dx)\\ &=O(n^{3\over 4}). \end{aligned} \]

若我们可以预处理一部分 \(S(k)\)，其中 \(k=1,2,\cdots,m\)，\(m\ge\lfloor\sqrt n\rfloor\)。设预处理时间复杂度为 \(T_0(m)\)，此时 \(T(n)\) 为：

\[\begin{aligned} T(n)&=T_0(m)+\sum_{k\in R(n),k>m}T(k)\\ &=T_0{m}+\sum_{k=1}^{\lfloor {n\over m}\rfloor}O(\sqrt{n\over k})\\ &=O(T_0(m)+\int_0^{n\over m}\sqrt{n\over x}\text dx)\\ &=O(T_0(m)+{n\over\sqrt m}). \end{aligned} \]

若使用线性筛预处理前缀和，则 \(T_0(m)=O(m)\)。由均值得：当 \(m=\Theta(n^{2\over3})\) 时，\(T(n)\) 取得最小值 \(O(n^{2\over 3})\)。

最后，对于多组询问要使用记忆化，防止重复求 \(O(\sqrt n)\) 种取值中的某一种。这样就可以保证复杂度是 \(O(n^{2\over3})\)。

应用

\(O(n)\) 跑不下来可以考虑尝试杜教筛。

筛 \(\mu\)

我们有恒等式

\[\mu \otimes1=\epsilon \]

启发我们取 \(g(n)=1(n)\)，\(h(n)=\epsilon(n)\)。所以有：

\[1(1)S(n)=\sum_{i=1}^n\epsilon(i)-\sum_{i=2}^n1(i)S(\lfloor{n\over i}\rfloor) \]

即：

\[S(n)=1-\sum_{i=2}^n1(i)S(\lfloor{n\over i}\rfloor) \]

\(1(i)\) 的前缀和即为 \(i\)，直接杜教筛复杂度 \(O(n^{2\over3})\)。

筛 \(\varphi\)

我们有恒等式

\[\varphi\otimes1=\text{id} \]

启发我们取 \(g(n)=1(n)\)，\(h(n)=\text{id}(n)\)。所以有：

\[1(1)S(n)=\sum_{i=1}^n\text{id}(i)-\sum_{i=2}^n1(i)S(\lfloor{n\over i}\rfloor) \]

即：

\[S(n)={n(n+1)\over2}-\sum_{i=2}^n1(i)S(\lfloor{n\over i}\rfloor) \]

\(1(i)\) 的前缀和即为 \(i\)，直接杜教筛复杂度 \(O(n^{2\over3})\)。

代码实现

Luogu P4213 【模板】杜教筛 Link

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

const int maxn = (1 << 21) + 10, maxv = (1 << 21);
ll T, n; 

vector<ll> primes, isp(maxn, true);
ll mu[maxn], phi[maxn];
unordered_map<int, ll> sphi;
unordered_map<int, int> smu;
void init() {
	mu[1] = 1, phi[1] = 1;
	for(int i = 2; i <= maxv; i++) {
		if(isp[i]) {
			primes.push_back(i);
			mu[i] = -1, phi[i] = i - 1;
		}
		for(int p : primes) {
			if(i * p > maxv) break;
			isp[i * p] = false;
			if(i % p) {
				mu[i * p] = -mu[i], phi[i * p] = phi[i] * (p - 1);
			}
			else {
				mu[i * p] = 0, phi[i * p] = phi[i] * p;
				break;
			}
		}
	}
	for(int i = 1; i <= maxv; i++) mu[i] += mu[i - 1], phi[i] += phi[i - 1];
	return;
}

ll get_phi(int x) {
	if(x <= maxv) return phi[x];
	if(sphi[x]) return sphi[x];
	ll res = 1ll * x * (1ll * x + 1) / 2;
	for(ll l = 2, r = 0; l <= x; l = r + 1) {
		r = x / (x / l);
		res -= get_phi(x / l) * (r - l + 1);
	}
	return sphi[x] = res;
}
int get_mu(int x) {
	if(x <= maxv) return mu[x];
	if(smu[x]) return smu[x];
	ll res = 1;
	for(ll l = 2, r = 0; l <= x; l = r + 1) {
		r = x / (x / l);
		res -= get_mu(x / l) * (r - l + 1);
	}
	return smu[x] = res;
}

int main() {
	ios :: sync_with_stdio(false); cin.tie(0); cout.tie(0);
	init();
	cin >> T; 
	while(T--) {
		cin >> n;
		cout << get_phi(n) << " " << get_mu(n) << endl;
	}
	return 0;
}

完结撒花 😃