积性函数的线性求法
首先我们来了解一下积性函数的定义。
\(\forall x,y\)满足\(gcd(x,y)=1\),有\(f(xy)=f(x)f(y)\)的数论函数称为积性函数。若\(\forall x,y\)均有\(f(xy)=f(x)f(y)\),则称为完全积性函数。
然后我们要理解线性素数筛。
void getP(int n)
{
for(int i=2;i<=n;i++)
{
if(!notP[i]) P[++cntP]=i;
for(int j=1;j<=cntP&&i*P[j]<=n;j++)
{
notP[i*P[j]]=true;
if(i%P[j]==0) break;
}
}
for(int i=1;i<=cntP;i++) printf("%d ",P[i]);
}
if(i%P[j]==0) break;
保证了每个数只会被其最小的质因子筛掉一次,所以复杂度为\(O(n)\)。用反证法证明:假设\(x\)最小的质因子为\(p_0\),而\(x\)是被\(i\times p_j(p_j>p_0)\)筛掉的,那么必然有\(p_0|i\),在循环到\(j\)之前已经break
了,矛盾。
那么接下来我们就可以来看线性求积性函数的代码了,这里以欧拉函数\(\varphi(n)\)为例。
void getF(int n)
{
f[1]=1,p1[1]=1;
for(int i=2;i<=n;i++)
{
if(!notP[i])
{
P[++cntP]=i;
for(lint j=i;j<=n;j*=i) f[j]=j-j/i,p1[j]=j; //1
}
for(int j=1;j<=cntP&&i*P[j]<=n;j++)
{
int x=i*P[j]; notP[x]=true;
if(i%P[j]) f[x]=f[i]*f[P[j]],p1[x]=P[j];
else {f[x]=f[i/p1[i]]*f[p1[i]*P[j]],p1[x]=p1[i]*P[j]; break;}
}
}
for(int i=1;i<=n;i++) printf("%lld ",f[i]); puts("");
}
将正整数\(x\)表示为\(\prod_{i=1}^k p_i^{q_i}\)的形式,其中\(p_i\)是递增的。
\(cntP\)与\(P[N]\)存储质数,\(f[x]\)记录函数的值,\(p1[x]\)记录\(x\)的\(p_1^{q_1}\)。
主要分析13,14行。
if(i%P[j]) f[x]=f[i]*f[P[j]],p1[x]=P[j];
因为\(i\geq P_j\)且\(i \ mod \ P_j \neq 0\),所以\(gcd(i,P_j)=1\),\(f(i\cdot P_j)=f(i)f(P_j)\);\(P_j\)是\(i\cdot P_j\)最小的质因子且只有一个,所以\(p1[i\cdot P_j]=P_j\)。
else {f[x]=f[i/p1[i]]*f[p1[i]*P[j]],p1[x]=p1[i]*P[j]; break;}
易知\(f(x\cdot p_1)=f(\dfrac{x}{p_1^{q_1}})f(p_1^{q_1+1})\),而\(P_j\)恰是\(i\)的\(p_1\);\(i\cdot P_j=p_1 \prod_{t=1}^k p_t^{q_t}=p_1^{q_1+1}\prod_{i=2}^k p_i^{q_i}\),所以\(p1[i\cdot P_j]=p_1^{q_1+1}=p1[i]\cdot P_j\)。本行的break
和素数筛一样保证了复杂度为\(O(n)\)。
不过当\(x=p_1^{q_1}\)时无法使用这个公式,所以要将注释1处的代码替换掉。