知识点 - 杜教筛

知识点 - 杜教筛

解决问题类型：

积性函数前缀和，

做法是 将原函数和恒等函数卷积得到F，对F求前缀和并使用微扰法求出原函数前缀和。

模板：

//注意mp<int,ll>
//一般很卡常，可以考虑换语言交
ll GetSum(int n) { // 算 f 前缀和的函数
  ll ans = f_g_sum(n); // 算 f * g 的前缀和
  // 以下这个 for 循环是整除分块
  for(ll l = 2, r; l <= n; l = r + 1) { // 注意从 2 开始
    r = (n / (n / l)); 
    ans -= (g_sum(r) - g_sum(l - 1)) * GetSum(n / l);
    // g_sum 是 g 的前缀和
    // 递归 GetSum 求解
  } return ans; 
}

公式推导：

问题:设f(n)为数论函数，求
$$S(n)=\sum _{i=1}^{n}f(n)$$
想法：根据$f(n)$的性质，构造一个$S(n)$关于$S(\lfloor \frac{n}{i}\rfloor )$的递推式。

构造方法如下：

找到一个合适的$g(n),\mathrm{let}F=f\ast g$ 即
$$F(n)=\sum _{d|n}f(\frac{n}{d})g(d)$$
对F求前缀和，则有
$$\begin{gathered} \sum F(n)=\sum _{i=1}^n\sum _{d|i}f(\frac{i}{d})g(d)=\sum _{d|i}\sum _{d=1}^{n}f(\frac{i}{d})g(d) \\ =\sum _{d=1}^{n}g(d)\sum _{i=1}^{\lfloor \frac{n}{d}\rfloor }f(i)=\sum _{d=1}^{n}g(d)S(\lfloor \frac{n}{d}\rfloor ) \end{gathered}$$
即
$$\sum _{i=1}^n\sum _{d|i}f(d)g(\frac{i}{d})=\sum _{i=1}^{n}g(i)S(\lfloor \frac{n}{i}\rfloor )$$
提出$g(1)S(n)$,移项（即微扰法 perturb the sum）可以得到递推式
$$g(1)S(n)=\sum _{i=1}^{n}F(i)-\sum _{i=2}^{n}g(i)S(\lfloor \frac{n}{i}\rfloor )$$
注意到后半部分的式子中的$S(\lfloor \frac{n}{i}\rfloor )$可以用整除分块$O(\sqrt{n})$地求出（见知识点 - 整除分块）。

复杂度

若能够快速求：（具体来说，小于$O(\sqrt{n})$）
$$\sum _{i=1}^n(f\ast g)(i)，\sum _{i=1}^{n}g(i)$$
那么计算单个S(n)的时间复杂度就变为
$$T(n)=O(\sqrt{n})+\sum _{i=1}^{\sqrt{n}}(T(i)+T(\frac{n}{i}))$$
展开一阶小量，积分，有
$$T(n)=\sum _{i=1}^{\sqrt{n}}O(\sqrt{i})+O(\sqrt{\frac{n}{i}})=O(n^{\frac{3}{4}})$$
可以打表预处理来加速，如果线性预处理出了f(n)的前k项$(k\ge \sqrt{n})$，此时复杂度变为：
$$T(n)=k+\sum _{i=1}^{\frac{n}{k}}\sqrt{\frac{n}{i}}=O(\frac{n}{\sqrt{k}}+k)$$
$k=n^{\frac{2}{3}}$时取到最小值$O(2n^{\frac{2}{3}})$

例题1

令$\Phi (n)={\sum }_{i=1}^n\phi (i)$为欧拉函数前缀和。

解析

利用性质
$$\sum _{d|n}\phi (d)=n\iff I\ast \phi =id$$
所以令$g(n)$为常函数$I(n)=1$,带入通式
$$g(1)S(n)=\sum _{i=1}^{n}F(i)-\sum _{i=2}^{n}g(i)S(\lfloor \frac{n}{i}\rfloor )$$
得到：
$$1\cdot \Phi (n)=\sum _{i=1}^{n}i-\sum _{i=2}^{n}1\cdot \Phi (\lfloor \frac{n}{i}\rfloor )$$
（补充推导：当然也可以从$I\ast \phi =id$直接从头推导，有
$$\therefore \phi (n)=n-\sum _{d|n,d<n}\phi (d)$$
推导得：
$$\begin{gathered} \varphi (n)=\sum _{i=1}^n\phi (i)=\sum _{i=1}^{n}i-\sum _{d|i,d<i}\phi (d)=\frac{n\cdot (n+1)}{2}-\sum _{i=2}^n\sum _{d|i,d<i}\phi (d)=\frac{n\cdot (n+1)}{2}-\sum _{\frac{i}{d}=2}^n\sum _{d=1}^{\lfloor \frac{n}{\frac{i}{d}}\rfloor }\phi (d) \\ =\frac{n\cdot (n+1)}{2}-\sum _{i=2}^n\sum _{d=1}^{\lfloor \frac{n}{i}\rfloor }\phi (d)=\frac{n\cdot (n+1)}{2}-\sum _{i=2}^n\varphi (\lfloor \frac{n}{i}\rfloor ) \end{gathered}$$
）

理论分析

具体数学4.9节的做法：

首先对于实数x有结论：
$$\begin{gathered} \sum _{d\ge 1}\Phi (\frac{x}{d})=\lfloor x\rfloor \cdot \lfloor x+1\rfloor /2 \\ \mathrm{where}\Phi (x)=\Phi (\lfloor x\rfloor ) \end{gathered}$$
（微扰法将第1项提出来有
$$\Phi (x)=\lfloor x\rfloor \cdot \lfloor x+1\rfloor /2-\sum _{d\ge 2}\Phi (\lfloor \frac{x}{d}\rfloor )$$
​ 得到杜教筛的递推式。）

第一个式子可以这么理解：可以将右边看成所有$0\le m<n\le x$的分数$\frac{m}{n}$的（未约分）分数的数量，而左边可以理解为形式为$\frac{m}{n}$且$\gcd (m,n)=d$的分数的个数为$\Phi (\frac{x}{d})$.

精确推导如下：
$$\frac{n\cdot (n+1)}{2}=\sum _{i=1}^{n}i=\sum _{i=1}^n\sum _d[d|i]\cdot \phi (d)=\sum _{\frac{i}{d}=1}^n\sum _{d=1}^{\lfloor \frac{n}{\frac{i}{d}}\rfloor }\phi (d)=\sum _{i=1}^n\varphi (\lfloor \frac{n}{i}\rfloor )$$
直接对第一个式子使用莫比乌斯反演变式:
$$g(x)=\sum _{d\ge 1}f(x/d)\iff f(x)=\sum _{d\ge 1}\mu (d)g(x/d)$$
注意这个变式与原莫比乌斯公式的区别：
$$g(n)=\sum _{d|n}f(d)\iff f(n)=\sum _{d|n}g(\frac{n}{d})\mu (d)$$
该变式适用于所有
$$\sum _{k,d\ge 1}|f(x/kd)|<\infty$$
这里顺便证明一下：
$$\begin{gathered} \sum _{d\ge 1}\mu (d)g(x/d)=\sum _{d\ge 1}\mu (d)\sum _{k\ge 1}f(x/kd) \\ =\sum _{m\ge 1}f(x/m)\sum _{k,d\ge 1}\mu (d)[m=kd] \\ =\sum _{m\ge 1}f(x/m)\sum _{d|m}\mu (d)=\sum _{m\ge 1}f(x/m)[m=1]=f(x) \end{gathered}$$
${\sum }_{d\ge 1}\Phi (\frac{x}{d})$用反演，得到
$$\Phi (x)=1/2\sum _{d\ge 1}\mu (d)\lfloor x/d\rfloor \cdot \lfloor x/d+1\rfloor$$
由于$\lfloor x/d\rfloor$前有莫比乌斯函数的系数，无法用整除分块，复杂度依然为O(n)。

不过这个式子可以用来计算近似值:(具体数学ch9中有证明)
$$\Phi (x)=3/{\pi }^2{x}^2+O(xlogx)$$

例题2

求梅滕斯函数$M(n)={\sum }_{i=1}^n\mu (i)$

解析

利用性质
$$[n=1]=\sum _{d|n}\mu (d)n\iff I\ast \mu =e$$
和上一题同样地，令$g(n)$为常函数$I(n)=1$,带入通式
$$g(1)S(n)=\sum _{i=1}^{n}F(i)-\sum _{i=2}^{n}g(i)S(\lfloor \frac{n}{i}\rfloor )$$
得到
$$1\cdot \mathrm{M}(n)=1-\sum _{i=2}^{n}1\cdot \mathrm{M}(\lfloor \frac{n}{i}\rfloor )$$
（补充推导：
$$1=\sum _{i=1}^n[i=1]=\sum _{i=1}^n\sum _{d|i}\mu (d)=\sum _{i=1}^n\sum _{d=1}^{\lfloor \frac{n}{i}\rfloor }\mu (d)=\sum _{i=1}^{n}M(\lfloor \frac{n}{i}\rfloor )$$
）

代码：

P4213 【模板】杜教筛（Sum）

#include<iostream>
#include<map>
#include<string>
#include<cstring>
#include<vector>
#include<algorithm>
#include<set>
#include<sstream>
#include<cstdio>
#include<cmath>
#include<climits>
#include<cstdlib>
#include<iomanip>
#include<unordered_map>
using namespace std;
const double PI = acos(-1.0);

#define rep(i,j,k) for(int i = (int)j;i <= (int)k;i ++)
#define debug(x) cerr<<#x<<":"<<x<<endl
#define pb push_back
#define FAST_IO ios::sync_with_stdio(0); cin.tie(0)
const int maxn = 1e7;
const int N = 3e6;
typedef long long ll;
int mu[N], lp[N], sumu[N];
ll phi[N],sump[N];
vector<int> pr;
void init() {
	mu[1] = 1;
	phi[1] = 1;
	for (int i = 2; i <= N; ++i) {
		if (lp[i] == 0) {
			lp[i] = i;
			mu[i] = -1;
			phi[i] = i - 1;
			pr.push_back(i);
		}
		for (int j = 0; j < (int)pr.size() && pr[j] <= lp[i] && i*pr[j] <= N; ++j)
		{
			lp[i * pr[j]] = pr[j];
			if (i%pr[j])phi[i*pr[j]] = phi[i] * (pr[j] - 1), mu[i * pr[j]] = -mu[i];
			else phi[i*pr[j]] = phi[i] * pr[j];
		}
	}
	sumu[0] = 0;
	sump[0] = 0;
	rep(i, 1, N-1) {
		sumu[i] = sumu[i-1] + mu[i];
		sump[i] = sump[i - 1] + phi[i];
	}
}
unordered_map<int,ll> Smu, Sphi;
inline ll GetSumu(int n) {
	if (n <= N-1) return sumu[n]; // sumu是提前筛好的前缀和
	if (Smu[n]) return Smu[n]; // 记忆化
	ll ret = 1ll; // 单位元的前缀和就是 1
	for (int l = 2, r; l <= n; l = r + 1) {
		r = n / (n / l); ret -= (r - l + 1) * GetSumu(n / l);
		// (r - l + 1) 就是 I 在 [l, r] 的和
	} return Smu[n] = ret; // 记忆化
}
inline ll GetSphi(int n) {
	if (n <= N-1) return sump[n]; // 提前筛好的
	if (Sphi[n]) return Sphi[n]; // 记忆化
	ll ret = 1ll * n * (n + 1) / 2; // f * g = id 的前缀和
	for (int l = 2, r; l <= n; l = r + 1) {
		r = n / (n / l); ret -= (r - l + 1) * GetSphi(n / l);
		// 同上，因为两个的 g 都是 I 
	} return Sphi[n] = ret; // 记忆化
}
int main() {
	init();
	int t;
	cin >> t;
	while (t--) {
		int n;
		cin >> n;
		cout << GetSphi(n) << ' ' << GetSumu(n) << endl;
	}
	int n;
	cin >> n;
}
/*
6
1
2
8
13
30
2333
*/