素數算法補充之"篩法"

時間 2020-10-02

標籤素數算法補充篩法简体版

原文原文鏈接

國慶中秋雙節，就不寫太長的文章了。node

補充和複習一下之前沒寫的素數區間篩法算法吧算法

部分證實過程來自OI wiki數組

素數篩法

若是咱們想要知道小於等於 \(n\) 有多少個素數呢？函數

一個天然的想法是咱們對於小於等於 \(n\) 的每一個數進行一次斷定。這種暴力的作法顯然不能達到最優複雜度，考慮如何優化。優化

考慮這樣一件事情：若是 \(x\) 是合數，那麼 \(x\) 的倍數也必定是合數。利用這個結論，咱們能夠避免不少次沒必要要的檢測。spa

若是咱們從小到大考慮每一個數，而後同時把當前這個數的全部（比本身大的）倍數記爲合數，那麼運行結束的時候沒有被標記的數就是素數了。code

int Eratosthenes(int n) {
    int p = 0;
    for (int i = 0; i <= n; ++i) is_prime[i] = 1;
    is_prime[0] = is_prime[1] = 0;
    for (int i = 2; i <= n; ++i) {
        if (is_prime[i]) {
        prime[p++] = i;  // prime[p]是i,後置自增運算表明當前素數數量
        // 由於從 2 到 i - 1 的倍數咱們以前篩過了，這裏直接從 i
        // 的倍數開始，提升了運行速度
        for (int j = i * i; j <= n;j += i)  is_prime[j] = 0;  
    }
    return p;
}

以上爲 Eratosthenes 篩法 （埃拉託斯特尼篩法），時間複雜度是 \(O(n\log\log n)\) 。it

以上作法仍有優化空間，咱們發現這裏面彷佛會對某些數標記了不少次其爲合數。有沒有什麼辦法省掉無心義的步驟呢？table

答案固然是：有！class

若是能讓每一個合數都只被標記一次，那麼時間複雜度就能夠降到 \(O(n)\) 了

void init() {
    phi[1] = 1;
    for (int i = 2; i < MAXN; ++i) {
        if (!vis[i]) {
            phi[i] = i - 1;
            pri[cnt++] = i;
        }
        for (int j = 0; j < cnt; ++j) {
            if (1ll * i * pri[j] >= MAXN) break;
            vis[i * pri[j]] = 1;
            if (i % pri[j]) {
                phi[i * pri[j]] = phi[i] * (pri[j] - 1);
            }
            else {
                // i % pri[j] == 0
                // 換言之，i 以前被 pri[j] 篩過了
                // 因爲 pri 裏面質數是從小到大的，因此 i 乘上其餘的質數的結果必定也是
                // pri[j] 的倍數 它們都被篩過了，就不須要再篩了，因此這裏直接 break
                // 掉就行了
                phi[i * pri[j]] = phi[i] * pri[j];
                break;
            }
        }
    }
}

上面代碼中的 \(phi\) 數組，會在下面提到。

上面的這種 線性篩法 也稱爲 Euler 篩法 （歐拉篩法）。node: 注意到篩法求素數的同時也獲得了每一個數的最小質因子

篩法求歐拉函數

注意到在線性篩中，每個合數都是被最小的質因子篩掉。好比設 \(p_1\) 是 \(n\) 的最小質因子， \(n' = \frac{n}{p_1}\) ，那麼線性篩的過程當中 \(n\) 經過 \(n' \times p_1\) 篩掉。

觀察線性篩的過程，咱們還須要處理兩個部分，下面對 \(n' \bmod p_1\) 分狀況討論。

若是 \(n' \bmod p_1 = 0\) ，那麼 \(n'\) 包含了 \(n\) 的全部質因子。

\[\begin{aligned} \varphi(n) & = n \times \prod_{i = 1}^s{\frac{p_i - 1}{p_i}} \\\\ & = p_1 \times n' \times \prod_{i = 1}^s{\frac{p_i - 1}{p_i}} \\\\ & = p_1 \times \varphi(n') \end{aligned} \]

那若是 \(n' \bmod p_1 \neq 0\) 呢，這時 \(n'\) 和 \(p_1\) 是互質的，根據歐拉函數性質，咱們有：

\[\begin{aligned} \varphi(n) & = \varphi(p_1) \times \varphi(n') \\\\ & = (p_1 - 1) \times \varphi(n') \end{aligned} \]

void phi_table(int n, int* phi) {
    for (int i = 2; i <= n; i++) phi[i] = 0;
    phi[1] = 1;
    for (int i = 2; i <= n; i++)
        if (!phi[i])
            for (int j = i; j <= n; j += i) {
                if (!phi[j]) phi[j] = j;
                phi[j] = phi[j] / i * (i - 1);
            }
}

篩法求莫比烏斯函數

線性篩

void pre() {
    mu[1] = 1;
    for (int i = 2; i <= 1e7; ++i) {
        if (!v[i]) mu[i] = -1, p[++tot] = i;
        for (int j = 1; j <= tot && i <= 1e7 / p[j]; ++j) {
            v[i * p[j]] = 1;
            if (i % p[j] == 0) {
                mu[i * p[j]] = 0;
                break;
            }
            mu[i * p[j]] = -mu[i];
        }
    }

篩法求約數個數

用 \(d_i\) 表示 \(i\) 的約數個數， \(num_i\) 表示 \(i\) 的最小質因子出現次數。

約數個數定理

定理：若 \(n=\prod_{i=1}^mp_i^{c_i}\) 則 \(d_i=\prod_{i=1}^mc_i+1\) .

證實：咱們知道 \(p_i^{c_i}\) 的約數有 \(p_i^0,p_i^1,\dots ,p_i^{c_i}\) 共 \(c_i+1\) 個，根據乘法原理， \(n\) 的約數個數就是 \(\prod_{i=1}^mc_i+1\) .

實現

由於 \(d_i\) 是積性函數，因此可使用線性篩。

void pre() {
    d[1] = 1;
    for (int i = 2; i <= n; ++i) {
        if (!v[i]) v[i] = 1, p[++tot] = i, d[i] = 2, num[i] = 1;
        for (int j = 1; j <= tot && i <= n / p[j]; ++j) {
            v[p[j] * i] = 1;
            if (i % p[j] == 0) {
                num[i * p[j]] = num[i] + 1;
                d[i * p[j]] = d[i] / num[i * p[j]] * (num[i * p[j]] + 1);
                break;
            }
            else {
                num[i * p[j]] = 1;
                d[i * p[j]] = d[i] * 2;
            }
        }
    }
}

篩法求約數和

\(f_i\) 表示 \(i\) 的約數和， \(g_i\) 表示 \(i\) 的最小質因子的 \(p+p^1+p^2+\dots p^k\) .

void pre() {
    g[1] = f[1] = 1;
    for (int i = 2; i <= n; ++i) {
        if (!v[i]) v[i] = 1, p[++tot] = i, g[i] = i + 1, f[i] = i + 1;
        for (int j = 1; j <= tot && i <= n / p[j]; ++j) {
            v[p[j] * i] = 1;
            if (i % p[j] == 0) {
                g[i * p[j]] = g[i] * p[j] + 1;
                f[i * p[j]] = f[i] / g[i] * g[i * p[j]];
                break;
            }
            else {
                f[i * p[j]] = f[i] * f[p[j]];
                g[i * p[j]] = 1 + p[j];
            }
        }
    }
    for (int i = 1; i <= n; ++i) f[i] = (f[i - 1] + f[i]) % Mod;
}

區間篩法

給定整數a和b，請問區間\([a,b)\)內有多少個素數？\(（a<b≤10^{12},b-a≤10^6）\)

思路：b之內的合數的最小質因數必定不超過\(√b\)。若是有\(√b\)之內的素數表的話，就能夠把埃氏篩法運用在\([a,b)\)上了。也就是說，先分別作好\([2,√b)\)的表和\([a,b)\)的表，而後從\([2,√b)\)的表中篩得素數的同時，也將其倍數從\([a,b)\)的表中劃去，最後剩下的就是\([a,b)\)內的素數了。

bool v1[Max_n1]; //數組大小爲sqrt(b)
bool v2[Max_n2]; //數組大小爲b-a

ll Prime(ll a, ll b) {
    for (ll i = 0; i * i < b; i++)v1[i] = true;
    for (ll i = 0; i < b - a; i++)v2[i] = true;

    for (ll i = 2; i * i < b; i++) {
        if (v1[i]) {
            for (ll j = 2 * i; j * j < b; j += i)v1[j] = false; //篩[2,b)
            for (ll j = max(2LL, (a + i - 1) / i) * i; j < b; j += i)v2[j - a] = false; //篩[a,b)
            //2LL是2的長整數形式
            //((a+i-1)/i)*i是符合>=a最小是i倍數的數
        }
    }
    ll k = 0;
    for (ll i = 0; i < b - a; i++) {
        if (v2[i])k++;
    }
    return k;