积性函数的线性求法

首先我们来了解一下积性函数的定义。

\(\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处的代码替换掉。

posted @ 2018-04-07 23:36  VisJiao  阅读(310)  评论(0编辑  收藏  举报