欧拉筛求积性函数的一般形式

【普通的欧拉筛】

int prime[MAXN/10], cntprime, f[MAXN];
bool vis[MAXN];
inline void sieve() {
	f[1]=1;
	for(int i=2; i<=Lim; ++i) {
		if(!vis[i]) {
			prime[++cntprime]=i;
			f[i]=...;
		}
		for(int j=1; j<=cntprime; ++j)
			if(prime[j]*i>Lim) break;
			else {
				vis[ prime[j]*i ]=1;
				if(i%prime[j]) 
					f[i*prime[j]]=f[i]*f[prime[j]];
				else {
					f[i*prime[j]]=...;
					break;
				}
			}
	}
}

由于没有维护最小质因子,一般用于线筛一些比较简单的积性函数


【维护最小质因子的欧拉筛】

\(fc_n\) 表示 \(n\) 的最小质因子

int fc[MAXN], prime[MAXN/10], cntprime, f[MAXN];
inline void sieve() {
	f[1]=1;
	for(int i=2; i<=Lim; ++i) {
		if(!fc[i]) {
			fc[i]=prime[++cntprime]=i;
			f[i]=...;
		}
		for(int j=1; j<=cntprime; ++j)
			if(prime[j]*i>Lim) break;
			else {
				fc[ prime[j]*i ]=prime[j];
				if( fc[i]!=prime[j] )
					f[ prime[j]*i ]=...;
				else {
					f[ prime[j]*i ]=...;
					break;
				}
			}
	}
}

由于维护了最小质因数,方便线筛一些 \(\boldsymbol f(p^{k+1})\)\(\boldsymbol f(p^k)\) 有较好递推关系的积性函数

且由于将运算过程中,用于判定的最小质因数保存下来,避免重复取模,一定程度上可以提高运行效率

但和前者一样,对于例如除数和函数 \(\sigma(n)\) 等,在取模意义下时,由于递推关系需要用到乘除法,可能会因为快速幂等退化为 \(O(n\log P)\)


【维护最小质因子次数的欧拉筛】

\(fcv_n\) 表示 \(n\) 的所有最小质因子乘积,\(fck_n\) 表示 \(n\) 的最小质因子次数(有的函数需要用到,没用到的可以删去)

int fcv[MAXN], fck[MAXN], prime[MAXN/10], cntprime, f[MAXN];
inline void sieve() {
	f[1]=1;
	for(int i=2; i<=Lim; ++i) {
		if(!fcv[i]) {
			prime[++cntprime]=i;
			fcv[i]=i;
			fck[i]=1;
			f[i]=...;
		}
		for(int j=1; j<=cntprime; ++j)
			if(prime[j]*i>Lim) break;
			else if(fcv[i]%prime[j]) {
				fcv[ i*prime[j] ]=prime[j];
				fck[ i*prime[j] ]=1;
				f[ i*prime[j] ]=f[i]*f[ prime[j] ];
			}
			else if(fcv[i]==i) {
				fcv[ i*prime[j] ]=prime[j]*i;
				fck[ i*prime[j] ]=fck[i]+1;
				f[ i*prime[j] ]=...;
				break;
			}
			else {
				fcv[ i*prime[j] ]=prime[j]*fcv[i];
				fck[ i*prime[j] ]=fck[i]+1;
				f[ i*prime[j] ]=f[i/fcv[i]]*f[ fcv[i]*prime[j] ];
				break;
			}
	}
}

由于运算过程中,每次转移等价于目标数去除最小质因子的结果,乘上最小质因子的贡献。避免了因为取模意义下需要求逆元消除贡献,而退化为 \(O(n\log P)\)

理论上该算法能在线性时间内处理比较复杂的积性函数,但维护的数据较多,且乘除法次数也较多,可能效率会相对降低


【维护最小质因子次数的欧拉筛-空间换“时间”】

在前者的基础上,采用空间换时间的操作,避免多次的乘除法:

维护 \(fc_n\) 表示 \(n\) 的最小质因数,避免检验最小质因数时采用取模检验

维护 \(lst_n\) 表示 \(n\) 去除最小质因数后的剩余部分,避免暴力除去最小质因数

int fcv[MAXN], fck[MAXN], prime[MAXN/10], cntprime, f[MAXN], fc[MAXN], lst[MAXN];
inline void sieve() {
	f[1]=1;
	for(int i=2; i<=Lim; ++i) {
		if(!fcv[i]) {
			prime[++cntprime]=i;
			fcv[i]=i;
			fck[i]=1;
			fc[i]=i;
			lst[i]=1;
			f[i]=...;
		}
		for(int j=1; j<=cntprime; ++j)
			if(prime[j]*i>Lim) break;
			else if(fc[i]>prime[j]) {
				fcv[ i*prime[j] ]=prime[j];
				fck[ i*prime[j] ]=1;
				fc[ i*prime[j] ]=prime[j];
				lst[ i*prime[j] ]=i;
				f[ i*prime[j] ]=f[i]*f[ prime[j] ];
			}
			else if(fcv[i]==i) {
				fcv[ i*prime[j] ]=prime[j]*i;
				fck[ i*prime[j] ]=fck[i]+1;
				fc[ i*prime[j] ]=prime[j];
				lst[ i*prime[j] ]=lst[i];
				f[ i*prime[j] ]=...;
				break;
			}
			else {
				fcv[ i*prime[j] ]=prime[j]*fcv[i];
				fck[ i*prime[j] ]=fck[i]+1;
				fc[ i*prime[j] ]=prime[j];
				lst[ i*prime[j] ]=lst[i];
				f[ i*prime[j] ]=f[ lst[i] ]*f[ fcv[i]*prime[j] ];
				break;
			}
	}
}

空间换时间可以避免过多的乘除法。但由于维护的信息过多,可能反而不如前者优秀


【例题】

求解 \(\displaystyle \sum_{i=1}^n \boldsymbol {\sigma_k}(i)\)

由于 \(\boldsymbol {\sigma_k}(n)=\sum_{d\mid n}d^k\)

\(\boldsymbol {\sigma_k}(p^{t+1})=\boldsymbol {\sigma_k}(p^t)\cdot p^k+1\)

可以仅对每个质数预处理 \(p^k\) ,复杂度为 \(O(\pi(n)\cdot \log P)\approx O({n\ln n}\cdot \log P)\approx O(n)\)

#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
typedef pair<ll, ll> pii;
typedef double db;
#define fi first
#define se second
const int Lim=1e7, MAXN=Lim+10, P=1e9+7;
inline int fpow(ll a, ll x) { ll ans=1; for(;x;x>>=1,a=a*a%P) if(x&1) ans=ans*a%P; return ans; }
int k, fcv[MAXN], prime[MAXN/10], cntprime, sigmak[MAXN], pk[MAXN/10];
inline void sieve(int lim=Lim) {
	sigmak[1]=1;
	for(int i=2; i<=lim; ++i) {
		if(!fcv[i]) {
			prime[++cntprime]=i;
			fcv[i]=i;
			sigmak[i]=((pk[cntprime]=fpow(i, k))+1)%P;
		}
		for(int j=1; j<=cntprime; ++j)
			if(prime[j]*i>lim) break;
			else if(fcv[i]%prime[j]) {
				fcv[ i*prime[j] ]=prime[j];
				sigmak[ i*prime[j] ]=(ll)sigmak[i]*sigmak[ prime[j] ]%P;
			}
			else if(fcv[i]==i) {
				fcv[ i*prime[j] ]=prime[j]*i;
				sigmak[ i*prime[j] ]=((ll)sigmak[i]*pk[j]+1)%P;
				break;
			}
			else {
				fcv[ i*prime[j] ]=prime[j]*fcv[i];
				sigmak[ i*prime[j] ]=(ll)sigmak[i/fcv[i]]*sigmak[ fcv[i]*prime[j] ]%P;
				break;
			}
	}
}
int main(){
    ios::sync_with_stdio(0);
    cin.tie(0); cout.tie(0);
	int n;
	ll res=0;
	cin>>n>>k;
	sieve(n);
	for(int i=1; i<=n; ++i) res+=sigmak[i];
	cout<<res%P;
    cout.flush();
    return 0;
}

目前 榜6

posted @ 2021-09-30 16:04  JustinRochester  阅读(135)  评论(0编辑  收藏  举报