积性函数的线性求法

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

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

转载于:https://www.cnblogs.com/VisJiao/p/8735603.html