[洛谷P4213]【模板】杜教筛（Sum）(杜教筛讲解)

题面

https://www.luogu.com.cn/problem/P4213

题解

前置知识

杜教筛

对于一些数论函数，杜教筛可以在低于线性的时间复杂度内求出其前缀和。

适用范围：对于数论函数f,g,h，如果g、h的前缀和能够\(O(1)\)求解，并且有

\[h=f{\times}g$$（其中${\times}$表示狄利克雷卷积） 那么，我们可以在低于线性的时间复杂度内求出f的前缀和。原理如下：（以下sf,sg,sh分别表示f,g,h的前缀和） $$h(n)={\sum\limits_{d1d2=n}}f(d1)g(d2)\]

\[{\sum\limits_{i=1}^{n}}h(i)={\sum\limits_{d1d2{\leq}n}}f(d1)g(d2) \]

\[sh(i)={\sum\limits_{d2=1}^{n}}g(d2){\sum\limits_{d1=1}^{{\lfloor}{\frac{n}{d1}}{\rfloor}}}f(d1) \]

\[={\sum\limits_{d=1}^{n}}g(d)sf({\lfloor}{\frac{n}{d}}{\rfloor}) \]

将d=1一项提出，得

\[g(1)sf(n)=sh(n)-{\sum\limits_{i=2}^{n}g(i)sf({\lfloor}{\frac{n}{i}}{\rfloor})} \]

为了方便，不妨设\(g(1)=1\)。那么有

\[sf(n)=sh(n)-{\sum\limits_{i=2}^{n}}g(i)sf({\lfloor}{\frac{n}{i}}{\rfloor}) \]

即可递归计算了。使用数论分块可以将计算\(sf(n)\)的时间降到\(O({\sqrt{n}})\)（不算递归下去的时间）；并且，递归计算的所有子问题sf(n')，一定可以写成\(n'={\lfloor}{\frac{n}{t}}{\rfloor}\)的形式，其中t是某个正整数。所以根据数论分块，记忆化后所有需要计算sf的下标只有\({\lfloor}{\frac{n}{i}}{\rfloor}(1{\leq}i{\leq}n)\)去重后的那些数，计算的总时间是它们的平方根之和，也就是

\[{\sum\limits_{i=1}^{{\sqrt{n}}}}{\sqrt{i}}+{\sum\limits_{i=1}^{\sqrt{n}}}{\sqrt{\frac{n}{i}}} \]

两边分别做积分近似，可以得到这样的时间复杂度为\(O(n^{\frac{3}{4}})\)。

如果函数f是积性函数，那么我们还可以进行进一步的优化：设置阈值M，将1~M的f值线性筛出来，如果计算sf(n)时\(n{\leq}M\)就直接调用这些值，否则再进行记忆化搜索。

考虑M的最优值。如果\(M<{\sqrt{n}}\)，总时间为

\[O(M)+{\sum\limits_{i=M+1}^{\sqrt{n}}}+{\sum\limits_{i=1}^{\sqrt{n}}}{\sqrt{\frac{n}{i}}} \]

发现第三项还是\(O(n^{\frac{3}{4}})\)，所以这样不行。

如果\(M{\geq}{\sqrt{n}}\)，总时间为

\[O(M)+{\sum\limits_{i=1}^{\frac{n}{M}}}{\sqrt{\frac{n}{i}}} \]

积分近似得第二项是\(O({\frac{n}{\sqrt{m}}})\)，所以置\(M=n^{\frac{2}{3}}\)的时候，总时间复杂度最小，为\(O(n^{\frac{2}{3}})\)。

代码

#include<bits/stdc++.h>

using namespace std;

#define N 2147483647
#define M 5000000
#define T 1300  
#define ll long long
#define rg register

inline ll read(){
	ll s = 0,ww = 1;
	char ch = getchar();
	while(ch < '0' || ch > '9'){if(ch == '-')ww = -1;ch = getchar();}
	while('0' <= ch && ch <= '9'){s = 10 * s + ch - '0';ch = getchar();}
	return s * ww;
}

inline void write(ll x){
	if(x < 0)putchar('-'),x = -x;
	if(x > 9)write(x / 10);
	putchar('0' + x % 10);
}

ll pri[M/5+5];
ll pn;
bool isp[M+5];
ll smu[M+5],sphi[M+5];

inline void Eular(){
	smu[1] = sphi[1] = 1;
	for(rg ll i = 2;i <= M;i++)isp[i] = 1;
	for(rg ll i = 2;i <= M;i++){
		if(isp[i])pri[++pn] = i,smu[i] = -1,sphi[i] = i - 1;
		for(rg ll j = 1;i * pri[j] <= M;j++){
			isp[i*pri[j]] = 0;
			if(i % pri[j])smu[i*pri[j]] = -smu[i],sphi[i*pri[j]] = sphi[i] * sphi[pri[j]];
			else{
				smu[i*pri[j]] = 0;
				sphi[i*pri[j]] = sphi[i] * pri[j];
				break;
			}
		}
	}
	for(rg ll i = 2;i <= M;i++)smu[i] += smu[i-1],sphi[i] += sphi[i-1]; 
}

ll n;

unordered_map<ll,ll>Hmu,Hphi;

inline ll sumphi(ll x){
	if(x <= M)return sphi[x]; //小于M的部分，预处理出来，使用数组存储 
	if(Hphi.count(x))return Hphi[x]; //大于M的部分，记忆化搜索，使用Hash表存储 
	ll ans = 1ll * (ll)x * ((ll)x + 1) >> 1;
	ll L,R = 1;
	while(R < x){
		L = R + 1,R = x / (x / L);
		ans -= 1ll * (R - L + 1) * sumphi(x / L);
	}
	return Hphi[x] = ans;
}

inline ll summu(ll x){
	if(x <= M)return smu[x]; //同上
	if(Hmu.count(x))return Hmu[x]; //同上
	ll ans = 1;
	ll L,R = 1;
	while(R < x){
		L = R + 1,R = x / (x / L);
		ans -= 1ll * (R - L + 1) * summu(x / L);
	}
	return Hmu[x] = ans;
}

int main(){
	Eular();
	ll c = read();
	while(c--){
		n = read();
		write(sumphi(n)),putchar(' '),write(summu(n)),putchar('\n');
	}	
	return 0;
}
posted @ 2020-02-10 17:23  coder66  阅读(196)  评论(0)  编辑  收藏  举报