数学筛法
有的时候怕忘记,写篇小博客记录一下。
线性筛素数
inline void init(int n)
{
for(int i=2;i<=n;i++)
{
if(!vis[i]) prime[++cnt]=i;
for(int j=1;j<=cnt && i*prime[j];j++)
{
vis[i*prime[j]]=1;
if(!(i%prime[j])) break;
}
}
return;
}
线性筛求欧拉函数
inline void init(int n)
{
memset(vis,0,sizeof(vis));
phi[1]=1;
for(int i=2;i<=n;i++)
{
if(!vis[i]) {prime[++cnt]=vis[i]=i,phi[i]=i-1;}
for(int j=1;j<=cnt && i*prime[j]<=n;j++)
{
vis[i*prime[j]]=prime[j];
phi[i*prime[j]]=phi[i]*(i%prime[j]?prime[j]-1:prime[j]);
}
}
for(int i=1;i<=n;i++) sum[i]=sum[i-1]+phi[i];
return;
}
线性筛求莫比乌斯函数
inline void init(int n)
{
mu[1]=1;
for(int i=2;i<=n;i++)
{
if(!vis[i]) {mu[i]=-1;prime[++cnt]=i;}
for(int j=1;j<=cnt && i*prime[j]<=n;j++)
{
vis[i*prime[j]]=1;
if(!(i%prime[j])) {mu[i*prime[j]]=0;break;}
else mu[i*prime[j]]=-mu[i];
}
}
for(int i=1;i<=n;i++) sum[i]=sum[i-1]+mu[i];
return;
}
线性筛求约数个数
用 \(d_i\) 表示 \(i\) 的约数个数,\(num_i\) 表示 \(i\) 的最小质因子出现次数。
inline void init(int n)
{
d[1]=1;
for(int i=2;i<=n;i++)
{
if(!vis[i])
{
prime[++cnt]=i;
d[i]=2,num[i]=1;
}
for(int j=1;j<=cnt && i*prime[j]<=n;j++)
{
vis[i*prime[j]]=1;
if(!(i%prime[j]))
{
num[i*prime[j]]=num[i]+1;
d[i*prime[j]]=d[i]/num[i*prime[j]]*(num[i*prime[j]]+1);
break;
}
else num[i*prime[j]]=1,d[i*prime[j]]=d[i]*2;
}
}
return;
}
线性筛求约数和
用 \(s_i\) 表示 \(i\) 的约数和,\(g_i\) 表示 \(i\) 的最小质因子的 \(p^0+p^1+p^2+\dots +p^k\).
inline void init(int n)
{
s[1]=g[1]=1;
for(int i=2;i<=n;i++)
{
if(!vis[i]) vis[i]=1,prime[++cnt]=i,g[i]=i+1,s[i]=i+1;
for(int j=1;j<=cnt && i*prime[j]<=n;j++)
{
vis[i*prime[j]]=1;
if(!(i%prime[j]))
{
g[i*prime[j]]=g[i]*prime[j]+1;
s[i*prime[j]]=s[i]/g[i]*g[i*prime[j]];
break;
}
else
{
s[i*prime[j]]=s[i]*s[prime[j]];
g[i*prime[j]]=prime[j]+1;
}
}
}
return;
}
杜教筛
对于数论函数 \(f\),杜教筛可以在低于线性时间的复杂度内计算 \(S(n)=\sum\limits_{i=1}^{n}f(i)\)。
原理
我们要想办法构造一个 \(S(n)\) 关于 \(S(\lfloor\frac{n}{i}\rfloor)\) 的递推式。
对于任意一个数论函数 \(g\),一定满足:
\[\begin{aligned}
\sum\limits_{i=1}^n(f*g)(i)&=\sum\limits_{i=1}^n\sum\limits_{d|i}g(d)f\left(\frac{i}{d}\right)\\
&=\sum\limits_{i=1}^n g(i)S\left(\left\lfloor\frac{n}{i}\right\rfloor\right)
\end{aligned}
\]
那么可以得到递推式:
\[\begin{aligned}
g(1)S(n)&=\sum\limits_{i=1}^n g(i)S\left(\left\lfloor\frac{n}{i}\right\rfloor\right)-\sum\limits_{i=2}^n g(i)S\left(\left\lfloor\frac{n}{i}\right\rfloor\right)\\
&=\sum\limits_{i=1}^n(f*g)(i)-\sum\limits_{i=2}^n g(i)S\left(\left\lfloor\frac{n}{i}\right\rfloor\right)
\end{aligned}
\]
假如我们可以构造恰当的数论函数 \(g\) 使得:
- 可以快速计算 \(\sum_{i=1}^n(f * g)(i)\)。
- 可以快速计算 g 的单点值,以用数论分块求解 \(\sum_{i=2}^ng(i)S\left(\left\lfloor\frac{n}{i}\right\rfloor\right)\)。
则我们可以在较短时间内求得 \(g(1)S(n)\)。
那么就得到 \(S(n)=\frac{\sum\limits_{i=1}^n(f*g)(i)-\sum\limits_{i=2}^n g(i)S\left(\left\lfloor\frac{n}{i}\right\rfloor\right)}{g(1)}\)。
时间复杂度
反正是 \(\Theta(n^{\frac{2}{3}})\)
例题1 P4213【模板】杜教筛(Sum)
点击查看代码
#include<bits/stdc++.h>
#define int long long
using namespace std;
const int MAXN=5e6+5;
int T,n,cnt,a[MAXN];
int prime[MAXN],vis[MAXN],mu[MAXN],phi[MAXN],summu[MAXN],sumphi[MAXN];
int f[MAXN],tot[MAXN];
map<int,int>ansmu;
map<int,int>ansphi;
inline void init(int n)
{
mu[1]=1; phi[1]=1;
for(int i=2;i<=n;i++)
{
if(!vis[i]) prime[++cnt]=i,mu[i]=-1,phi[i]=i-1;
for(int j=1;j<=cnt && i*prime[j]<=n;j++)
{
vis[i*prime[j]]=1;
if(!(i%prime[j])) {mu[i*prime[j]]=0;phi[i*prime[j]]=phi[i]*prime[j];break;}
else mu[i*prime[j]]=-mu[i],phi[i*prime[j]]=phi[i]*(prime[j]-1);
}
}
for(int i=1;i<=n;i++) summu[i]=summu[i-1]+mu[i];
for(int i=1;i<=n;i++) sumphi[i]=sumphi[i-1]+phi[i];
return;
}
inline int getmu(int x)
{
if(x<=MAXN) return summu[x];
if(ansmu[x]) return ansmu[x];
int ans=1;
for(int l=2,r;l<=x;l=r+1)
{
r=x/(x/l);
ans-=(r-l+1)*getmu(x/l);
}
return ansmu[x]=ans;
}
inline int getphi(int x)
{
if(x<=MAXN) return sumphi[x];
if(ansphi[x]) return ansphi[x];
int ans=(x*(x+1))/2ll;
for(int l=2,r;l<=x;l=r+1)
{
r=x/(x/l);
ans-=(r-l+1)*getphi(x/l);
}
return ansphi[x]=ans;
}
signed main()
{
ios_base::sync_with_stdio(false);
cin.tie(0),cout.tie(0);
cin>>T; init(MAXN);
while(T--)
{
cin>>n;
printf("%lld %lld\n",getphi(n),getmu(n));
}
return 0;
}