数论——质数筛法

最新推荐文章于 2023-11-01 17:11:55 发布

爱寂寞的时光

最新推荐文章于 2023-11-01 17:11:55 发布

阅读量122

点赞数

分类专栏： # 数论文章标签：算法线性代数

本文链接：https://blog.csdn.net/jiahonghao2002/article/details/113702557

版权

数论专栏收录该内容

15 篇文章 3 订阅

订阅专栏

质数筛法

朴素筛法

枚举素数即可，判断一个数 $x$ 是不是素数，从 $2$ 开始到 $\sqrt{x}$ 枚举他的约数即可。

因为一个合数 $x$ 可以写成是两个约数的乘积 $m n = x$ ，假设 $\gt \sqrt{x}$ ，那么 $\lt \sqrt{x}$ ，因此，对于任意的合数 $x$ ，必定存在一个约数 $\leq \sqrt{x}$ 。

class Solution {
public:
    bool isPrime(int x) {
        for (int i = 2; i * i <= x; ++i) {
            if (x % i == 0) {
                return false;
            }
        }
        return true;
    }

    int countPrimes(int n) {
        int ans = 0;
        for (int i = 2; i < n; ++i) {
            ans += isPrime(i);
        }
        return ans;
    }
};

时间复杂度为 $O(n\sqrt{n})$ ，空间复杂度为 $O (1)$ 。

埃氏筛

我们考虑这样一个事实：如果 $x$ 是质数，那么大于 $x$ 的 $x$ 的倍数 $2x,3x,\ldots$ 一定不是质数，因此我们可以从这里入手。又因为唯一素约数分解定理，一个合数必定可以写成一个小于他的素数和一个整数的乘积，因此我们只需要碰见一个素数，就关掉它的倍数即可。

class Solution {
public:
    int countPrimes(int n) {
        vector<int> isPrime(n, 1);
        int ans = 0;
        for (int i = 2; i < n; ++i) {
            if (isPrime[i]) {
                ans += 1;
                if ((long long)i * i < n) {
                    for (int j = i * i; j < n; j += i) {
                        isPrime[j] = 0;
                    }
                }
            }
        }
        return ans;
    }
};

假设 $x$ 是质数，应该从 $x * x$ 开始关掉，因为 $2x,3x,4x,5x,\ldots,x(x-1)$ ，已经被 $2,3,4,5,\ldots,x-1$ 关掉，不需要再考虑了。

开根优化

同时可以加上开根优化，对于枚举小于等于 $n$ 的素数，我们只需要筛到 $\sqrt{n}$ 即可，因为，如果有大于 $\sqrt{n}$ 的合数 $m$ ，他必定有小于 $\sqrt{m}$ 的素因子，因为 $\sqrt{m} \leq \sqrt{n}$ ，因此合数 $m$ 一定能够被关掉。

int primes[100000];

void eratosthenes(int n)
{
	int i = 0;
	for (; i <= n; i++)
		primes[i] = 1;
	primes[0] = primes[1] = 0;

	for (int i = 2; i * i <= n; i++)
	{
		if (primes[i])
		{
			int j = i;
			for (; i * j <= n; j++)
				primes[i * j] = 0;
		}
	}
}

奇数优化

除了偶素数2之外的偶数都是合数，因此跳过偶数只处理奇数。

因为素数一定是奇数，我们优化一下埃氏筛：

只在奇数范围内寻找素数
奇数乘以偶数还是偶数，因此没有必要标记，因为在第一步就已经排除了，所以我们只标记奇合数就好了。

var countPrimes = function(n) {
    var isCom= new Int8Array(n), b = Math.sqrt(n), r = n > 2 ? 1 : 0, t
    for(var i = 3; i < n; i += 2) {
        if(!isCom[i]) {
            r++
            if (i <= b) for(var j = i; (t = i * j) < n; j += 2) isCom[t] = 1
        }
    }
    return r
};

可以证明，埃氏筛法的时间复杂度为 $\log \log n)$ 。

线性筛

在上面的埃氏筛中，一个合数可能被关闭多次，我们就有了优化的方向，让一个合数只被关闭一次即可。

根据唯一素因子分解定理，一个合数等于一个素约数乘以另外一个约数，那么这么说，我们可以建立一个线性表，里面储存了我们到目前为止遇到的素数，然后迭代途中，遇到一个数，就依次关掉这个线性表中的元素与这个数的乘积即可。

可以证明，线性表primes中的元素一定是线性递增的。我们可以证明这个算法的正确性：

如果遇到的是素数，那么只能被他自己和1分解，而在这之前没有理由关掉这个素数，因此这个数一定是素数。
如果遇到一个合数，那么这个合数一定可以写作是，一个小于这个数的质数 $p$ ，和大于这个质数 $p$ 但小于这个数的一个数 $k$ 的乘积 $p k$ 。（有时候也可以是小于这个质数 $p$ 的 $k$ ，如果 $k$ 是质数，那么 $p$ 和 $k$ 互相交换一下就好，如果 $k$ 是合数，那么一定能在 $k$ 中找到一个素数，把它当做 $p$ 就好，总之，这两个数 $p$ 和 $k$ 一定能找到） $p$ 一定是在 $k$ 访问之前就已经被添加到素数列表中了，之后再访问 $k$ 的时候，就会关掉合数 $p k$ 了。

class Solution {
public:
    int countPrimes(int n) {
        vector<int> primes;
        vector<int> isPrime(n, 1);
        for (int i = 2; i < n; ++i) {
            if (isPrime[i]) {
                primes.push_back(i);
            }
            for (int j = 0; j < primes.size() && i * primes[j] < n; ++j) {
                isPrime[i * primes[j]] = 0;
            }
        }
        return primes.size();
    }
};

但是反过来想一想，这么做一个合数还是可能被访问多次啊，比如 $\times 2 = 5 \times 8$ ，当访问到20的时候，40会被关掉一次，当访问到8的时候，40又被关掉了一次，根本没有做到优化，能不能让一个合数被唯一一个素数确定呢？

答案是肯定的，就是最小素因子，即一个合数，一定存在一个最小素因子 $p$ ，和一个数 $k$ 的乘积 $p k$ ，并且在 $k$ 的分解中，没有比 $p$ 还小的素因子了。此时，一个合数就能被唯一一个素数所确定。

那我们怎么改进代码呢？添加一句话：

class Solution {
public:
    int countPrimes(int n) {
        vector<int> primes;
        vector<int> isPrime(n, 1);
        for (int i = 2; i < n; ++i) {
            if (isPrime[i]) {
                primes.push_back(i);
            }
            for (int j = 0; j < primes.size() && i * primes[j] < n; ++j) {
                isPrime[i * primes[j]] = 0;
                if (i % primes[j] == 0) {
                    break;
                }
            }
        }
        return primes.size();
    }
};

其中

if (i % primes[j] == 0) {
    break;
}

也就是说，如果遇到一个素数可以整除i了，那么就停止标记。这是为什么呢？假如继续标记 $i * p r i m e s [j + 1]$ 这个数，又因为 $i$ 可以被 $p r i m e s [j]$ 整除，也就是说 $i * p r i m e s [j + 1] = p r i m e s [j] * n * p r i m e s [j + 1]$ ，其中 $n = i / p r i m e s [j]$ ，并且 $n$ 不能再被分解为比 $p r i m e s [j]$ 更小的素数了（如果能分解，那么 $i$ 还能被比 $p r i m e s [j]$ 更小的质数整除了）因此，数 $p r i m e s [j] * n * p r i m e s [j + 1]$ 分解中 $p r i m e s [j]$ 就是最小的素因子了，而 $p r i m e s [j + 1]$ 不是，所以不能再这里标记，应该在数 $n * p r i m e s [j + 1]$ 的时候标记（这个数大于 $i$ ，肯定能在之后访问到）。其他数也同理。

总结一句话就是，针对被关闭的数 $k = i * p r i m e s [j]$ ，必须令 $p r i m e s [j]$ 是其最小的素因子。

时间复杂度 $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;
      }
    }
  }
}

对于： $i$ 是素数的情况，很显然 $\varphi(i)=i-1$ 。

对于：phi[i * pri[j]] = phi[i] * (pri[j] - 1);来说， $i$ 与 $p r i [j]$ 互质，因此利用积性函数的性质求解即可。

对于：phi[i * pri[j]] = phi[i] * pri[j];来说， $p r i [j]$ 是 $i$ 的一个因子，我们把数进行质因数分解 $i * pri[j] = pri[j]^{m}k$ ，因为两者互质，因此

$\varphi(pri[j]^{m}k) = \varphi(pri[j]^{m}) \varphi(k) = (pri[j]^{m} - pri[j]^{m-1}) \varphi(k) = pri[j] (pri[j]^{m-1} - pri[j]^{m-2}) \varphi(k) = pri[j] \varphi(i)$

因此该算法是正确的。

求莫比乌斯函数

根据莫比乌斯函数的定义：

当 $n = 1$ 的时候 $\mu(n)=1$ ；当 $n=p_{1}p_{2} \ldots p_{m}$ 的时候，即任意多个质数相乘并且没有幂次质数， $\mu(n)=(-1)^{m}$ ；其他情况 $\mu(n)=0$


#define MAXN 1000

int phi[MAXN];
bool vis[MAXN];
int pri[MAXN];
int mu[MAXN];
int cnt = 0;

void euler()
{
    phi[1] = 1;
    mu[1] = 1;
    for (int i = 2; i < MAXN; ++i)
    {
        if (!vis[i])
        {
            phi[i] = i - 1;
            // 唯一一个质数为-1
            mu[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);
                // 利用莫比乌斯函数的积性
                mu[i * pri[j]] = -mu[i];
            }
            else
            {
                // i % pri[j] == 0
                // 换言之，i 之前被 pri[j] 筛过了
                // 由于 pri 里面质数是从小到大的，所以 i 乘上其他的质数的结果一定也是
                // pri[j] 的倍数 它们都被筛过了，就不需要再筛了，所以这里直接 break
                // 掉就好了
                phi[i * pri[j]] = phi[i] * pri[j];
                // pri[j] 是 i的倍数，存在幂次质数，函数值为0
                mu[i * pri[j]] = 0;
                break;
            }
        }
    }
}

求约数个数

除法函数：

$\sigma_{k}(n)=\sum_{d \mid n}d^{k}$

为积性函数。

则当 $k = 0$ 的时候，除法函数退化为因数个数函数。仍为积性函数。

#define MAXN 1000

bool vis[MAXN];
int pri[MAXN];
int pnxt[MAXN]; // 表示去掉所有质因子之后，剩下的质因子的乘积
int sigma[MAXN];
int cnt = 0;

void euler()
{
    pnxt[1] = 1;
    sigma[1] = 1;
    for (int i = 2; i < MAXN; ++i)
    {
        if (!vis[i])
        {
            // 质数的约数只有自己和1
            sigma[i] = 2;
            // 质数的最小质因子的数量为1，去掉后乘积为1
            pnxt[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])
            {
                // 积性
                sigma[i * pri[j]] = sigma[i] * 2;
                // i * pri[j]  只有一个最小质数因子为 pri[j]
                pnxt[i * pri[j]] = i;
            }
            else
            {
                // 添加了一个最小质因子
                pnxt[i * pri[j]] = pnxt[i];
                // pri[j] 是 i 的一个素约数
                sigma[i * pri[j]] = sigma[i] + sigma[pnxt[i]];
                break;
            }
        }
    }
}

对于sigma[i * pri[j]] = sigma[i] + sigma[pnxt[i]];来说：

$\sigma(i * p) = \sigma(p^{m}k) = \sigma(p^{m}) \sigma(k) = (m + 1) \sigma(k) = (\sigma(p^{m-1}) + 1)) \sigma(k) = \sigma(i) + \sigma(k)$

其中 $k$ 为去掉所有最小质因子后剩下的乘积。

求约数和

#define MAXN 1000

bool vis[MAXN];
int pri[MAXN];
int pnxt[MAXN]; // 表示去掉所有质因子之后，剩下的质因子的乘积
int psum[MAXN]; // 表示最下质因子的幂次和
int sigma1[MAXN];
int cnt = 0;

void euler()
{
    pnxt[1] = 1;
    sigma1[1] = 1;
    for (int i = 2; i < MAXN; ++i)
    {
        if (!vis[i])
        {
            // 质数的约数只有自己和1
            sigma1[i] = 1 + i;
            // 质数的最小质因子的数量为1，去掉后乘积为1
            pnxt[i] = 1;
            // p0 + p1
            psum[i] = 1 + i;
            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])
            {
                // 积性
                sigma1[i * pri[j]] = sigma1[i] * (1 + pri[j]);
                // i * pri[j]  只有一个最小质数因子为 pri[j]
                pnxt[i * pri[j]] = i;
                psum[i * pri[j]] = pri[j] + 1;
            }
            else
            {
                // 添加了一个最小质因子
                pnxt[i * pri[j]] = pnxt[i];
                // 将所有项都递推一步 并 加上 p0
                psum[i * pri[j]] = psum[i] * pri[j] + 1;
                // pri[j] 是 i 的一个素约数
                sigma1[i * pri[j]] = psum[i * pri[j]] * sigma1[pnxt[i]];
                break;
            }
        }
    }
}

对于sigma1[i * pri[j]] = psum[i * pri[j]] * sigma1[pnxt[i]];来说：

$\sigma_{1}(i * p) = \sigma_{1}(p^{m}) \sigma_{1}(k) = \sigma_{1}(k) \sum_{i = 0}^{m}p^{i}$

字母的含义同上。

求最大质因子


void euler()
{
    for (int i = 2; i < 5000005; i++)
    {
        if (!isNprime[i])
        {
            smf[i] = i;
            // 素数的最大质因子是其本身
            primes[tot++] = i;
        }

        for (int j = 0; j < tot && i * primes[j] < 5000005; j++)
        {
            int val = i * primes[j];
            isNprime[val] = true;
            // 我们知道，线性筛的质数是从小到大添加的，因此最大质因子一定是原来数i的最大质因子
            smf[val] = smf[i];
            if (i % primes[j] == 0)
            {
                break;
            }
        }
    }
}

求最小质因子

void euler()
{
    for (int i = 2; i < 5000005; i++)
    {
        if (!isNprime[i])
        {
            smf[i] = i;
            // 素数的最小质因子是其本身
            primes[tot++] = i;
        }

        for (int j = 0; j < tot && i * primes[j] < 5000005; j++)
        {
            int val = i * primes[j];
            isNprime[val] = true;
            // 我们知道，线性筛的质数是从小到大添加的，因此最小质因子一定是新添加的素数primes[j]
            smf[val] = primes[j];
            if (i % primes[j] == 0)
            {
                break;
            }
        }
    }
}

爱寂寞的时光

关注

0
点赞
踩
2

收藏

觉得还不错? 一键收藏
0
评论
数论——质数筛法

数论-质数筛法朴素筛法枚举素数即可，判断一个数xxx是不是素数，从222开始到x\sqrt{x}x枚举他的约数即可。因为一个合数xxx可以写成是两个约数的乘积mn=xmn=xmn=x，假设m>xm \gt \sqrt{x}m>x，那么n<xn \lt \sqrt{x}n<x，因此，对于任意的合数xxx，必定存在一个约数k≤xk \leq \sqrt{x}k≤x。class Solution {public: bool isPrime(int x) {
复制链接

扫一扫

专栏目录