bzoj3944: Sum（杜教筛）

题目

Description

求$S(n)={\sum }_{i=1}^{n}f(i)$

Solution

part1（核心公式推导）:

对于任何函数$g(i)$，都有：
${\sum }_{i=1}^n(f\ast g)(i)$
$={\sum }_{i=1}^n{\sum }_{d|i}f(d)g(\frac{i}{d})$
$={\sum }_{i=1}^{n}g(i){\sum }_{i\ast j\le n}f(j)$
$={\sum }_{i=1}^{n}g(i)S(\lfloor \frac{n}{i}\rfloor )$
那么$g(1)S(n)={\sum }_{i=1}^n(f\ast g)(i)-{\sum }_{i=2}^{n}g(i)S(\lfloor \frac{n}{i}\rfloor )$
杜教筛就要构造一个$g(i)$，使得${\sum }_{i=1}^n(f\ast g)(i)$与$g(i)$能快速计算

part2（复杂度证明）:

假设计算$S(n)$的复杂度为$T(n)$
那么$T(n)=O(\sqrt{n})+{\sum }_{i=1}^{\sqrt{n}}T(i)+T(\frac{n}{i})$（$i$从$1$枚举到$n$时所有不同的$T(\lfloor \frac{n}{i}\rfloor )$之和）
据$tangjz$大佬所说：这里只展开一层就可以了，更深层的复杂度是高阶小量
所以$T(n)={\sum }_{i=1}^{\sqrt{n}}O(\sqrt{i})+O(\sqrt{\frac{n}{i}})$
$={\sum }_{i=1}^{\sqrt{n}}O(\sqrt{\frac{n}{i}})$
$=\sqrt{n}{\int }_1^{\sqrt{n}}{x}^{-\frac{1}{2}}dx$
$=\sqrt{n}\cdot 2x^{\frac{1}{2}}{|}_1^{\sqrt{n}}$
$=O(n^{\frac{3}{4}})$
如果预处理了前$k(k\ge \sqrt{n})$个数的前缀和，那么
$T(n)=O(k)+{\sum }_{i=1}^{\frac{n}{k}}O(\sqrt{\frac{n}{i}})=O(k+\frac{n}{\sqrt{k}})$
当$k=O(n^{\frac{2}{3}})$时$T(n)=O(n^{\frac{2}{3}})$

part3（example）:

就拿$\phi (i)$举例吧，有一个公式$\phi \ast 1=i{d}_1$
所以$1(1)S(n)={\sum }_{i=1}^{n}i{d}_1(i)-{\sum }_{i=2}^{n}1(i)S(\lfloor \frac{n}{i}\rfloor )$
$S(n)=\frac{n(n+1)}{2}-{\sum }_{i=2}^{n}S(\lfloor \frac{n}{i}\rfloor )$

part4（代码实现中的注意点）:

$1.$因为数据不止$1$组，所以如果预处理$k$个数的话，复杂度为$O(k+T\frac{n}{\sqrt{k}})$
根据基本不等式，$k$取$(Tn{)}^{\frac{2}{3}}$最为合适
$2.$计算$\frac{n(n+1)}{2}$时，因为有$n+1$，所以可能会爆$int$，数论分块中的$l$和$r$也是如此
$3.$关于记忆化，有两种方法。一种是开unordered_map，另外一种是假设现在正在计算$S(n^{\prime })$，那么只要判断$p[\lfloor \frac{n}{n^{\prime }}\rfloor ]$即可，数组大小可开到$O(\frac{n}{k})的规模$

Code

#include<cstdio>
#include<cstring>
typedef unsigned long long ll;
typedef unsigned int ui;
const int N=5000000;
int T,i,pr[N/14],n,cnt,j,t,mu[N],m[1291],nn;
bool vis[N];
ll phi[N],p[1291];
ll calp(int n){
	if (n<N) return phi[n];
	if (p[nn/n]) return p[nn/n];
	ll ans=(1ll*n+1)*n>>1;
	for (ui l=2,r;l<=n;l=r+1) r=n/(n/l),ans-=calp(n/l)*(r-l+1);
	return p[nn/n]=ans;
}
int calm(int n){
	if (n<N) return mu[n];
	if (m[nn/n]) return m[nn/n];
	int ans=1;
	for (ui l=2,r;l<=n;l=r+1) r=n/(n/l),ans-=calm(n/l)*(r-l+1);
	return m[nn/n]=ans;
}
int main(){
	scanf("%d",&T);
	phi[1]=mu[1]=1;
	for (i=2;i<N;i++){
		if (!vis[i]) pr[cnt++]=i,phi[i]=i-1,mu[i]=-1;
		for (j=0;j<cnt && (t=i*pr[j])<N;j++){
			vis[t]=1;
			if (!(i%pr[j])){
				phi[t]=phi[i]*pr[j];
				mu[t]=0;
				break;
			}
			phi[t]=phi[i]*(pr[j]-1);
			mu[t]=-mu[i];
		}
	}
	for (i=2;i<N;i++) phi[i]+=phi[i-1],mu[i]+=mu[i-1];
	for (;T--;) scanf("%d",&n),nn=n,memset(p,0,sizeof(p)),
	memset(m,0,sizeof(m)),printf("%llu %d\n",calp(n),calm(n));
}