【学习笔记】杜教筛

如果我们要求一个积性函数 f(x) 的前缀和,可以用杜教筛在 O(n23) 的复杂度求出。
具体地,构造函数 g(x) 和函数 h(x) ,使得 h=fg,要求的式子是 S(n)=i=1nf(i)
开始推式子。

i=1nh(i)=i=1ng(i)i=1d|iig(d)×f(id)i=1nh(i)=i=1ng(i)j=1njf(j)i=1nh(i)=g(1)S(n)+=i=1ng(i)×S(ni)i=1nh(i)=g(1)S(n)+i=2ng(i)×S(ni)g(1)S(n)=i=1nh(i)i=2ng(i)×S(ni)

如果h(x)的前缀和好做,后面的式子可以整除分快。
 
对于莫比乌斯函数,可以构造 g(x)I函数。
那么:S(n)=1i=2S(nd)
 
对于欧拉函数,可以构造为 I.
那么: S(n)=i=1nid=2nS(nd)

 
对于S(n)=i=1ni×ϕ(i),可以构造函数 g(x)=nd

对于某一类积性函数,可以通过观察积性函数的关系来构造函数。

 

P4213 【模板】杜教筛(Sum)代码
#include<bits/stdc++.h>
using namespace std;
#define int long long
inline int rd(){
	int f=1,j=0;
	char w=getchar();
	while(!isdigit(w)){
		if(w=='-')f=-1;
		w=getchar();
	}
	while(isdigit(w)){
		j=j*10+w-'0';
		w=getchar();
	}
	return f*j;
}
const int N=2e7+10,ma=2e7;
unordered_map<int,long long>mul,fi;
long long Mul[N],Fi[N];
void work(int t){
	if(t>ma&&mul.count(t))return ;
	if(t<=ma)return ;
	long long ansa=1ll*t*(t+1)/2,ansb=1;
	for(int l=2,r;l<=t;l=r+1){
		r=t/(t/l),work(t/l);
		if(t/l>ma)ansa-=1ll*(r-l+1)*mul[t/l],ansb-=1ll*(r-l+1)*fi[t/l];
		else ansa-=1ll*(r-l+1)*Mul[t/l],ansb-=1ll*(r-l+1)*Fi[t/l];
	}
	mul[t]=ansa,fi[t]=ansb;
	return ;
}
bool wk[N];
int zhi[N],tail;
void init(){
	Mul[1]=Fi[1]=1;
	for(int i=2;i<=ma;i++){
		if(!wk[i])zhi[++tail]=i,Mul[i]=i-1,Fi[i]=-1;
		for(int j=1;j<=tail&&zhi[j]*i<=ma;j++){
			wk[zhi[j]*i]=true;
			if(i%zhi[j]==0){
				Mul[i*zhi[j]]=Mul[i]*zhi[j];
				break;
			}
			Mul[i*zhi[j]]=Mul[zhi[j]]*Mul[i];
			Fi[i*zhi[j]]=-Fi[i];
		}
	}
	for(int i=2;i<=ma;i++)Mul[i]+=Mul[i-1],Fi[i]+=Fi[i-1];
	return ;
}
signed main(){
	init();
	int T=rd();
	while(T--){
		int t=rd();
		work(t);
		if(t>ma)printf("%lld %lld\n",mul[t],fi[t]);
		else printf("%lld %lld\n",Mul[t],Fi[t]);
	}
	return 0;
}
有点小卡常,要先线性筛出前 2e7 位。
posted @   flywatre  阅读(10)  评论(0编辑  收藏  举报
相关博文:
阅读排行:
· 地球OL攻略 —— 某应届生求职总结
· 周边上新:园子的第一款马克杯温暖上架
· Open-Sora 2.0 重磅开源!
· 提示词工程——AI应用必不可少的技术
· .NET周刊【3月第1期 2025-03-02】
点击右上角即可分享
微信分享提示