[数论] 欧拉筛法

· · 个人记录

前言

最近学数论，我是真的绝望，欧拉筛法也只能靠背代码勉强凑合凑合，但在我社CSQ大佬的帮助下，我理解到了其中神奇的奥妙

正题

欧拉筛法是一种可以筛出质数，欧拉函数，约数个数和约数和的筛法 那么我们就对这些问题逐一进行讲解

在这之前，我们先说几个东西：

1、每一个大于等于2的正整数n，都有n=p_1^{w_1}p_2^{w_2}…p_m^{w_m}（p_1到p_m按升序排列）

2、正整数n的欧拉函数phi(n)=n(1-\frac{1}{p_1})(1-\frac{1}{p_2})…(1-\frac{1}{p_m})=p_1^{w_1-1}(p_1-1)*p_2^{w_2-1}(p_2-1)*…*p_m^{w_m-1}(p_m-1)

3、正整数n的约数个数d(n)=(1+w_1)(1+w_2)…(1+w_m)

4、正整数n的约数和s(n)=(1+p_1+p_1^2+…+p_1^{w_1})(1+p_2+p_2^2+…+p_2^{w_2})……(1+p_m+p_m^2+…+p_m^{w_m})

质数

思路

筛质数是一个相对简单的内容，在普通的筛法中，我们知道当一个数i为质数时，可以筛掉所有i的倍数

但在欧拉筛法中，我们会进行一定的优化，那么外层循环枚举i，判断i是否为质数；内层循环枚举j，这个j实质上就是枚举第j个质数，筛掉i*prim[j]，重点就是接下来的步骤：

如果i\%prim[j]==0，我们就要break掉，这是因为我们筛素数是肯定是希望每一个数都能被它的最小质因子筛掉，当i\%prim[j]==0时，就相当于i=prim[j]*k（k为一个正整数），那么如果我们继续往后枚举到了一个质数prim[m]，就要筛去i*prim[m]，也就是prim[j]*k*prim[m]，因为我们枚举的质数是从小到大的，所以i*prim[m]的最小质因子应该是prim[j]，若继续枚举，就会出现冗余的计算

代码

inline void sieve(int x) {
    for(reg int i = 2;i <= x;i ++) {
        if(! vis[i])
            prim[++ len] = i;
        for(reg int j = 1;j <= len && i * prim[j] <= x;j ++) {
            vis[i * prim[j]] = 1;
            if(i % prim[j] == 0)
                break;
        }
    }
}

欧拉函数

思路

我们知道欧拉函数是这样的：phi(n)=p_1^{w_1-1}(p_1-1)*p_2^{w_2-1}(p_2-1)*…*p_m^{w_m-1}(p_m-1)

那么在之前筛质数的基础上，我们知道会进行两种判断，我们也根据这两种情况来求解：

1、i\%prim[j]==0，这意味着i里面包含了至少一个prim[j]，所以i*prim[j]也包含了至少一个prim[j]，且i*prim[j]在i的基础上多包含了一个prim[j]

在i的质因数分解里，p_1就代表prim[j]

因为phi(i)=p_1^{w_1-1}(p_1-1)*p_2^{w_2-1}(p_2-1)*…*p_m^{w_m-1}(p_m-1) 所以phi(i*prim[j])=p_1^{w_1-1+1}(p_1-1)*p_2^{w_2-1}(p_2-1)*…*p_m^{w_m-1}(p_m-1)

可得phi(i*prim[j])=phi(i)*p_1，即phi(i*prim[j])=phi(i)*prim[j]

2、i\%prim[j]!=0，就是说i里面没有包含prim[j]，所以i*prim[j]里在i的基础上只包含了一个prim[j]

我们把prim[j]记作p_{m+1}，指数为1

因为phi(i)=p_1^{w_1-1}(p_1-1)*p_2^{w_2-1}(p_2-1)*…*p_m^{w_m-1}(p_m-1) 所以phi(i*prim[j])=p_1^{w_1-1}(p_1-1)*p_2^{w_2-1}(p_2-1)*…*p_m^{w_m-1}(p_m-1)*p_{m+1}^{1-1}(p_{m+1}-1)

可得phi(i*prim[j])=phi(i)*(p_{m+1}-1)，即phi(i*prim[j])=phi(i)*(prim[j]-1)

代码

inline void sieve(int x) {
    phi[1] = 1;
    for(reg int i = 2;i <= x;i ++) {
        if(! vis[i]) {
            prim[++ len] = i;
            phi[i] = i - 1; //因为欧拉函数代表小于这个数的且与这个数互质的数的个数，所以质数的欧拉函数为它本身减1
        }
        for(reg int j = 1;j <= len && i * prim[j] <= x;j ++) {
            vis[i * prim[j]] = 1;
            if(i % prim[j] == 0) {
                phi[i * prim[j]] = phi[i] * prim[j];
                break;
            }
            phi[i * prim[j]] = phi[i] * (prim[j] - 1);
        }
    }
}

约数个数

思路

先摆公式d(n)=(1+w_1)(1+w_2)…(1+w_m)

和之前的一样，我们分两种情况讨论：

1、i\%prim[j]==0，此时i里面包含了至少一个prim[j]，i*prim[j]里面也至少包含了一个prim[j]

在i的质因数分解中，就有p_1=prim[j]

因为d(i)=(1+w_1)(1+w_2)…(1+w_m) 所以d(i*prim[j])=(1+w_1+1)(1+w_2)…(1+w_m)

可得d(i*prim[j])=d(i)/(1+prim[j])*(2+prim[j])

但是因为我们在用欧拉筛法求解时，并不知道w_1的值，所以我们用一个sum数组来计算并保存这个值，sum[i]表示i的最小质因数的个数

2、i\%prim[j]!=0，此时i里面没有包含prim[j]这个质因子，所以i*prim[j]里面只包含了一个prim[j]

也就是说，i=p_1^{w_1}p_2^{w_2}……p_m^{w_m}，i*prim[j]=p_1^{w_1}p_2^{w_2}……p_m^{w_m}p_{m+1}^{w_{m+1}}，而且我们知道w_{m+1}=1

因为d(i)=(1+w_1)(1+w_2)…(1+w_k)…(1+w_m) 所以d(i*prim[j])=(1+w_1)(1+w_2)…(1+w_k)…(1+w_m)(1+w_{m+1})

可得d(i*prim[j])=d(i)*2

代码

inline void sieve(int x) {
    for(reg int i = 2;i <= x;i ++) {
        if(! vis[i]) {
            prim[++ len] = i;
            d[i] = 2;   //质数的约数只有1和它本身
            sum[i] = 1;
        }
        for(reg int j = 1;j <= len && i * prim[j] <= x;j ++) {
            vis[i * prim[j]] = 1;
            if(i % prim[j] == 0) {
                sum[i * prim[j]] = sum[i] + 1;
                d[i * prim[j]] = d[i] / (sum[i] + 1) * (sum[i] + 2);
                break;
            }
            sum[i * prim[j]] = 1;
            d[i * prim[j]] = d[i] * 2;
        }
    }
}

约数和

思路

公式在此：s(n)=(1+p_1+p_1^2+……+p_1^{w_1})(1+p_2+p_2^2+……+p_2^{w_2})……(1+p_m+p_m^2+……+p_m^{w_m})

同样，我们分两种情况：

1、i\%prim[j]==0，就是i里面至少包含了一个prim[j]，所以i*prim[j]里面也包含了至少一个prim[j]，而且在i*prim[j]里面的prim[j]比i多了一个

在i的分解中，prim[j]实质上就是p_1，则

i=p_1^{w_1}p_2^{w_2}…p_m^{w_m}$，$i*prim[j]=p_1^{w_1+1}p_2^{w_2}…p_k^{w_k+1}…p_m^{w_m}

因为s(i)=(1+p_1+p_1^2+……+p_1^{w_1})(1+p_2+p_2^2+……+p_2^{w_2})…(1+p_m+p_m^2+……+p_m^{w_m}) 所以s(i*prim[j])=(1+p_1+p_1^2+……+p_1^{w_1}+p_1^{w_1+1})(1+p_2+p_2^2+……+p_2^{w_2})…(1+p_m+p_m^2+……+p_m^{w_m})

此时我们用一个psum数组，其中psum[i]表示i的最小质因子p的一个形如1+p+p^2+……+p^w的等比数列的值 那么psum[i*prim[j]] = psum[i]*prim[j]+1，s(i*prim[j])=s(i)/psum[i]*psum[i*prim[j]]

2、i\%prim[j]!=0，此时i里面没有包含prim[j]，那么i*prim[j]里面就只包含了一个prim[j]

假设i=p_1^{w_1}p_2^{w_2}…p_m^{w_m}，那么i*prim[j]=p_1^{w_1}p_2^{w_2}…p_m^{w_m}p_{m+1}^{w_{m+1}}

因为s(i)=(1+p_1+p_1^2+……+p_1^{w_1})(1+p_2+p_2^2+……+p_2^{w_2})…(1+p_m+p_m^2+……+p_m^{w_m}) 所以s(i*prim[j])=(1+p_1+p_1^2+……+p_1^{w_1})(1+p_2+p_2^2+……+p_2^{w_2})…(1+p_m+p_m^2+……+p_m^{w_m})(1+p_{m+1})

可得s(i*prim[j])=s(i)*(1+prim[j])

代码

inline void sieve(int x) {
    for(reg int i = 2;i <= x;i ++) {
        if(! vis[i]) {
            prim[++ len] = i;
            psum[i] = s[i] = i + 1;
        }
        for(reg int j = 1;j <= len && i * prim[j] <= x;j ++) {
            vis[i * prim[j]] = 1;
            if(i % prim[j] == 0) {
                psum[i * prim[j]] = psum[i] * prim[j] + 1;
                s[i * prim[j]] = s[i] / psum[i] * psum[i * prim[j]]
                break;
            }
            psum[i * prim[j]] = prim[j] + 1;
            s[i * prim[j]] = s[i] * psum[i * prim[j]];
        }
    }
}