淺談積性函數的線性篩法

前置知識

數論函數及相關基本定義

素數的線性篩

線性篩

線性篩可以在嚴格$O(n)$的時間內篩出積性函數的值，

它有常見的套路

假設$n = p_1^{a_1} p_2^{a_2} \dots p_k^{a_k}$

如果我們能快速得出$f(p_i)$和$f(p_i^{k+1})$的取值，那么直接套板子就行了

在下文中如無特殊說明，默認$p_i$表示$n$質因數分解之后第$i$個質數，$a_i$表示$p_i$的指數

常見的有以下幾種

線性篩素數

比較簡單，這也是篩其他積性函數的基礎

#include<cstdio>
const int MAXN = 1e4 + 10;
int N, prime[MAXN], vis[MAXN], tot;
void GetPrime(int N) {
    vis[1] = 1;
    for(int i = 2; i <= N; i++) {
        if(!vis[i]) prime[++tot] = i;
        for(int j = 1; j <= tot && i * prime[j] <= N; j++) {
            vis[i * prime[j]] = 1;
            if(!(i % prime[j])) break;
        }
    }
}
int main() {
    GetPrime(1e3);
    for(int i = 1; i <= tot; i++)
        printf("%d ", prime[i]);
    return 0;
}

線性篩莫比烏斯函數

這個也是比較常見的

根據莫比烏斯函數的定義

$$\mu =\begin{cases}\left( -1\right) ^{k}\left( n=p_{1}p_{2}\ldots p_{k}\right) \\ 0\left( \exists P^{2}|n\right) \\ 1\left( n=1\right) \end{cases}$$

直接篩就可以了

#include<cstdio>
const int MAXN = 1e4 + 10;
int N, prime[MAXN], vis[MAXN], mu[MAXN], tot;
void GetMu(int N) {
    vis[1] = mu[1] = 1;
    for(int i = 2; i <= N; i++) {
        if(!vis[i]) prime[++tot] = i, mu[i] = -1;
        for(int j = 1; j <= tot && i * prime[j] <= N; j++) {
            vis[i * prime[j]] = 1;
            if(!(i % prime[j])) {
                mu[i * prime[j]] = 0;
                break;
            }
            mu[i * prime[j]] = mu[i] * mu[prime[j]];
            //根據莫比烏斯函數的定義，這里也可以寫為
            //mu[i * prime[j]] = -mu[i];
        }
    }
}
int main() {
    GetMu(1e3);
    for(int i = 1; i <= tot; i++)
        printf("%d ", mu[i]);
    return 0;
}

 

線性篩歐拉函數

這個貌似更常用一點qwq。

我在以前的文章中也詳細的講過，這里只放一下代碼

#include<cstdio>
const int MAXN = 1e4 + 10;
int N, prime[MAXN], vis[MAXN], phi[MAXN], tot;
void GetPhi(int N) {
    vis[1] = phi[1] = 1;
    for(int i = 2; i <= N; i++) {
        if(!vis[i]) prime[++tot] = i, phi[i] = i - 1;
        for(int j = 1; j <= tot && i * prime[j] <= N; j++) {
            vis[i * prime[j]] = 1;
            if(!(i % prime[j])) {
                phi[i * prime[j]] = phi[i] * prime[j];
                break;
            }
            phi[i * prime[j]] = phi[i] * phi[prime[j]];
        }
    }
}
int main() {
    GetPhi(1e3);
    for(int i = 1; i <= tot; i++)
        printf("%d ", phi[i]);
    return 0;
}

線性篩約數個數

這個就稍微高端一點了，按照上面的套路，我們只需要考慮最小的質因子對每個數的貢獻

$n = p_1^{a_1} p_2^{a_2} \dots p_k^{a_k}$

設$d(i)$表示$i$的約數個數

那么根據約數定理

$d(i) = \prod_{i = 1}^k (a_i + 1) $

記$a(i)$表示$n$的最小的質因子($a_1$)的指數。

$d(p_i) = 2$

當$i % p_j = 0$時，考慮$i * p_j$，實際上也就是$a_i$的指數多了$1$

我們先除去原來的，再加上新的就OK了

#include<cstdio>
const int MAXN = 1e4 + 10;
int N, prime[MAXN], vis[MAXN], D[MAXN], a[MAXN], tot;
void GetD(int N) {
    vis[1] = D[1] = a[1] = 1;
    for(int i = 2; i <= N; i++) {
        if(!vis[i]) prime[++tot] = i, D[i] = 2, a[i] = 1;
        for(int j = 1; j <= tot && i * prime[j] <= N; j++) {
            vis[i * prime[j]] = 1;
            if(!(i % prime[j])) {
                D[i * prime[j]] = D[i] / (a[i] + 1) * (a[i] + 2);
                a[i * prime[j]] = a[i] + 1;
                break;
            }
            D[i * prime[j]] = D[i] * D[prime[j]];
            a[i * prime[j]] = 1;
        }
    }
}
int main() {
    GetD(1e3);
    for(int i = 1; i <= tot; i++)
        printf("%d ", D[i]);
    return 0;
}

線性篩約數和

同樣，按照上面的套路考慮

$n = p_1^{a_1} p_2^{a_2} \dots p_k^{a_k}$

設$SD(i)$表示$i$的約數和

$SD(n) = \prod_{i = 1}^k (\sum_{j = 1}^{a_i} p_i^j)$

$sum(i)$表示$i$最小的質因子的貢獻，即$sum(i) = \sum_{i = 1}^{a_1}p_1^j$

$low(i)$表示$i$最小質因子的指數，$low(i) = a_1$

有了這三個我們就可以轉移了

同樣是考慮$i$的最小的因子對答案的貢獻，應該比較好推，看代碼吧

#include<cstdio>
const int MAXN = 1e4 + 10;
int N, prime[MAXN], vis[MAXN], SD[MAXN], sum[MAXN], low[MAXN], tot;
void GetSumD(int N) {
    vis[1] = SD[1] = low[1] = sum[1] = 1;
    for(int i = 2; i <= N; i++) {
        if(!vis[i]) prime[++tot] = i, sum[i] = SD[i] = i + 1, low[i] = i;
        for(int j = 1; j <= tot && i * prime[j] <= N; j++) {
            vis[i * prime[j]] = 1;
            if(!(i % prime[j])) {
                low[i * prime[j]] = low[i] * prime[j];
                sum[i * prime[j]] = sum[i] + low[i * prime[j]];
                SD[i * prime[j]] = SD[i] / sum[i] * sum[i * prime[j]];
                break;
            }
            low[i * prime[j]] = prime[j];
            sum[i * prime[j]] = prime[j] + 1;
            //這里low和sum不是積性函數 
            SD[i * prime[j]] = SD[i] * SD[prime[j]];
        }
    }
}
int main() {
    GetSumD(1e3);
    for(int i = 1; i <= tot; i++)
        printf("%d ", SD[i]);
    return 0;
}

 

非常見積性函數的篩法

很多情況下我們會遇到求兩個積性函數狄利克雷卷積的情況

很顯然，這個函數也是積性函數，我們考慮如何求得

為了方便篩，我們需要把問題無限簡化，

設$low(i)$表示$p_1^{a_1}$

考慮篩法中最關鍵的地方

$i \% p_j = 0$，

有了$low(i)$，此時我們需要分兩種情況討論

1. $low(i) = i$，此時$i$一定是某個素數的冪的形式（否則就會break掉）

這里就用到了我最開始說的那個套路

如果我們能快速的利用$f(p_i^{k})$更新出$f(p_i^{k + 1})$，那這個素數就很容易篩了

2. $low(i) \not = i$，那么$i / low(i)$一定與$low(i) * p_j$是互質的，我們可以直接利用積性函數的性質去更新

C++版的偽代碼

vis[1] = low[1] = 1; H[1] = 初始化 
for(int i = 2; i <= N; i++) {
    if(!vis[i]) prime[++tot] = i, mu[i] = -1, H[i] = 質數的情況, low[i] = i;
    for(int j = 1; j <= tot && i * prime[j] <= N; j++) {
        vis[i * prime[j]] = 1;
        if(!(i % prime[j])) {
            low[i * prime[j]] = (low[i] * prime[j]); 
            if(low[i] == i) H[i * prime[j]] = 特殊判斷;
            else H[i * prime[j]] = H[i / low[i]] * H[prime[j] * low[i]];
            break;
        } 
        H[i * prime[j]] = H[i] * H[prime[j]];
        low[i * prime[j]] = prime[j];
    }
}

 

參考資料

積性函數與線性篩

線性篩約數個數和、約數和

線性篩,積性函數,狄利克雷卷積,常見積性函數的篩法