2019.01.21【BZOJ3944】【洛谷P4213】Sum（杜教筛）

洛谷传送门

BZOJ传送门

---

解析：

为什么这道题不在zxyoi的模板分类那一栏呢？

因为zxyoi要更这道题的题解，但是他从不更新模板题解

因为数学题代码是没有模板的，只有其中的思路是相通的。

所以zxyoi这里就直接开始推了。

前置姿势：

整除分块，积性函数，$Dirichlet$卷积，莫比乌斯反演

杜教筛：

应该都知道杜教是谁吧，不知道的请自行百度。

所以杜教筛是个什么东西，由杜教发明的一种在低于线性时间内处理出积性函数前缀和的东西。

杜教筛要解决的问题就是求这样一个函数$S(n)$的值：$$S(n)=\sum _{i=1}^{n}f(i)$$

其中$f(i)$是积性函数。

考虑一般怎么求解，既然$f$是积性函数，那么我们就可以线性筛了啊，然后求出前缀和，不就$O(n)$了吗。

好的，当$n\le 1e11$呢？

我们这时候需要一个低于线性时间复杂度的方法，这就是杜教筛要做的事情。

解决：

我们现在找到另一个积性函数$g(i)$，把$f$和$g$做一下$Dirichlet$卷积，求前缀和：$$\begin{array}{rl}& \sum _{i=1}^n(f\ast g)(i)\\ =& \sum _{i=1}^n\sum _{d\mid i}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{array}$$

我们考虑求$g(1)S(n)$，通过基础差分可以知道：$$\sum _{i=1}^{n}g(d)S(\lfloor \frac{n}{d}\rfloor )-\sum _{i=2}^{n}g(d)S(\lfloor \frac{n}{d}\rfloor )$$

所以杜教筛基础就完了，我们得到一个递归式
$$g(1)S(n)=\sum _{i=1}^n(f\ast g)(i)-\sum _{i=2}^{n}g(d)S(\lfloor \frac{n}{d}\rfloor )$$

这个式子用整除分块+递归求解+记忆化来实现，就已经十分优秀了。

至于具体的复杂度是多少，我们来算一算：$$T(n)=\sum _{i=1}^{\sqrt{n}}O(\sqrt{i})+O(\sqrt{\frac{n}{i}})=O(n^{\frac{3}{4}})$$

这个可以把求和转换成积分后得到，详细的就不再叙述，以免偏离重点。

对于小部分的数据我们可以预处理来做到$O(n^{\frac{2}{3}})$，考虑呢预处理前$m$个之后：$$T(n)=\sum _{i=1}^{\lfloor \frac{n}{m}\rfloor }O(\sqrt{\frac{n}{i}})+O(m)=O(\sqrt{\frac{n}{m}}+O(m))$$

解出来得到$m=O(n^{\frac{2}{3}})$的时候最优，此时$T(n)=O(n^{\frac{2}{3}})$。

当然这是单个询问的做法，多组数据还是尽量能多预处理就多预处理吧

---

对于本题：

1.我们要求的是$S(n)={\sum }_{i=1}^n\varphi (n)$

$Dirichlet$卷积基础知识：$\varphi \ast 1=Id$

所以设$f=\varphi ,g=1,f\ast g=Id$

直接得到：$$S(n)=\frac{n(n+1)}{2}-\sum _{i=2}^{n}S(\lfloor \frac{n}{i}\rfloor )$$

直接做就行了

2.我们要求的是：$S(n)={\sum }_{i=1}^n\varphi (n)$

$Dirichlet$卷积基础知识：$\mu \ast 1=ϵ$

设$f=\mu ,g=1,f\ast g=ϵ$

得到：$$S(n)=1-\sum _{i=2}^{n}S(\lfloor \frac{n}{i}\rfloor )$$

---

然后讲几个这道题的坑点。

洛谷上面有点卡常。

BZOJ上面有$n=0$和$n=2^{31}-1$的数据，记得要处理一下。

---

代码（递归版，好写但常数略大）：

#include<bits/stdc++.h>
#include<tr1/unordered_map>
using namespace std;
#define ll long long
#define re register
#define gc getchar
#define pc putchar
#define cs const

tr1::unordered_map<int,int> Summu;
tr1::unordered_map<int,ll> Sumphi;
cs int P=8000007,N=P-7,INF=0x7fffffff;
int prime[P],pcnt;
bool mark[P]; 
int mu[P];
ll phi[P];
inline void linear_sieves(int len=N){
    phi[1]=mu[1]=1;
    for(int re i=2;i<=len;++i){
        if(!mark[i]){
            prime[++pcnt]=i;
            phi[i]=i-1;
            mu[i]=-1;
        }
        for(int re j=1;i*prime[j]<=len;++j){
            mark[i*prime[j]]=true;
            if(i%prime[j]==0){
                phi[i*prime[j]]=phi[i]*prime[j];
                mu[i*prime[j]]=0;
                break;
            }
            phi[i*prime[j]]=phi[i]*(prime[j]-1);
            mu[i*prime[j]]=-mu[i];
        }
        mu[i]+=mu[i-1];
        phi[i]+=phi[i-1];
    }
}

inline int Get_summu(int n){
    if(n<=N)return mu[n];
    if(Summu[n])return Summu[n];
    int ans=1;
    for(int re i=2,j;j<INF&&i<=n;i=j+1){
        j=n/(n/i);
        ans-=(j-i+1)*Get_summu(n/i);
    }
    return Summu[n]=ans;
}

inline ll Get_sumphi(int n){
    if(n<=N)return phi[n];
    if(Sumphi[n])return Sumphi[n];
    ll ans=n*((ll)n+1)/2;
    for(int re i=2,j;j<INF&&i<=n;i=j+1){
        j=n/(n/i);
        ans-=(j-i+1)*Get_sumphi(n/i);
    }
    return Sumphi[n]=ans;
}

int T;
int n;
signed main(){
    linear_sieves();
    for(scanf("%d",&T);T--;){
        scanf("%d",&n);
        cout<<Get_sumphi(n)<<" "<<Get_summu(n)<<"\n";
    }
    return 0;
}

代码（非递归版，对优点就是跑得贼快）：

#include<bits/stdc++.h>
using namespace std;
#define ll long long
#define re register
#define gc getchar
#define pc putchar
#define cs const

cs int P=1700006,N=P-6;
int prime[P],pcnt,mu[P];
bool mark[P];
ll phi[P];

inline void linear_sieves(int len=N){
	phi[1]=mu[1]=1;
	for(int re i=2;i<=N;++i){
		if(!mark[i])prime[++pcnt]=i,phi[i]=i-1,mu[i]=-1;
		for(int re j=1,k;(k=i*prime[j])<=len;++j){
			mark[k]=true;
			if(i%prime[j]==0){
				phi[k]=phi[i]*prime[j];
				mu[k]=0;
				break;
			}
			phi[k]=phi[i]*(prime[j]-1);
			mu[k]=-mu[i];
		}
		phi[i]+=phi[i-1];
		mu[i]+=mu[i-1]; 
	}
}

int G[P];
inline int Get_summu(int n){
	if(n==0)return 0;
	int v=n/(int)pow(n,2.0/3);
	int res,last,now;
	for(int re i=v;i;--i){
		int m=n/i;
		res=0;
		for(int re j=2;i*j<=v;++j)res+=G[i*j];
		re int j;
		for(j=v/i+1;(ll)j*j<=m;++j)res+=mu[m/j];
		last=j-1;
		j=m/j;
		for(;j;--j){
			now=m/j;
			res+=mu[j]*(now-last);
			last=now;
		}
		G[i]=1-res;
	}
	return G[1];
}

ll F[P];
inline ll Get_sumphi(int n){
	if(n==0)return 0;
	int v=n/(int)pow(n,2.0/3);
	ll res;int last,now;
	for(int re i=v;i;--i){
		int m=n/i;
		res=0;
		for(int re j=2;i*j<=v;++j)res+=F[i*j];
		int j;
		for(j=v/i+1;(ll)j*j<=m;++j)res+=phi[m/j];
		last=j-1;
		j=m/j;
		for(;j;--j){
			now=m/j;
			res+=phi[j]*(now-last);
			last=now;
		}
		F[i]=m*((ll)m+1)/2-res;
	}
	return F[1];
}

int T,n;
signed main(){
	linear_sieves();
	for(scanf("%d",&T);T--;){
		scanf("%d",&n);
		cout<<Get_sumphi(n)<<" "<<Get_summu(n)<<"\n";
	}
	return 0;
}