C语言丨筛法求素数（质数）

素数（质数）是指在大于1的自然数中，除了1和它本身以外不再有其他因数的自然数。素数被广泛用于密码学、汽车变速箱齿轮设计、害虫的生物生长周期与杀虫剂使用之间的关系、导弹和鱼雷等领域上，具有重要意义。本文就来介绍求素数的一种方法：筛法。

---

在初学编程时，我们解决问题的想法应该都是定义法。按照素数的定义——除了1和它本身以外不再有其他因数的大于1的自然数，我们可以这样判断一个整数x是否为素数：

首先，x(x>1)满足不能被1和x以外的其他数整除，即我们可以遍历2~x-1，看是否有整数能整除x，这种方法称为试商法。为此我们可以写出最初版的代码：

int IsPrime(int x)
{
    int i, flag = 1;
    if (x <= 1)//若x<=1，则不是素数
        flag = 0;
    for (i=2; i<x && flag; i++)
    {
        if (x%i == 0)//若x能被2~x-1中的整数整除，则不是素数
            flag = 0;
    }
    return flag;//若是素数则返回1，否则返回0
}

但是当x较大时，循环的次数较多，效率低下。实际上，只要在2~x-1中发现了任意一个整数能够整除x，则可立刻判断x不是素数，此时已经可以跳出循环，不再进行重复无意义的工作，即：

int IsPrime(int x)
{
    int i, flag = 1;
    if (x <= 1)//若x<=1，则不是素数
        flag = 0;
    for (i=2; i<x && flag; i++)
    {
        if (x%i == 0)//若x能被2~x-1中的整数整除，则不是素数
        {
            flag = 0;
            break;//只要找到第一个，即可跳出循环
        }
    }
    return flag;//若是素数则返回1，否则返回0
}

实际上，只需用2~sqrt(x)之间的整数去试商看是否能被整除即可（这里不给出数学证明，有兴趣的朋友可以自行查阅资料），为此我们可以把循环条件改为i<=sqrt(x)，此时若x不是素数，循环次数将大大减少，可以进一步改进我们的代码。

int IsPrime(int x)
{
    int i, flag = 1;
    int squareRoot = sqrt(x);
    if (x <= 1)//若x<=1，则不是素数
        flag = 0;
    for (i=2; i<=squareRoot && flag; i++)
    {
        if (x%i == 0)//若x能被2~sqrt(x)中的整数整除，则不是素数
        {
            flag = 0;
            break;//只要找到第一个，即可跳出循环
        }
    }
    return flag;//若是素数则返回1，否则返回0
}

若我们要实现求100以内的所有素数，可在主函数中调用IsPrime()函数如下：

int main()
{
    int i;
    for (i=1; i<=100; i++)
    {
        if (IsPrime(i))
            printf("%d\n", i);
    }
    return 0;
}

但是试想，如果我们要实现求素数的范围较大（如100000以内），那么在遍历1~100000的循环中，每一次循环都要调用一次IsPrime()函数，且在函数中还要遍历2~sqrt(i)，这样将会使程序的效率变得异常低下。事实上，我们完全可以找出较少的数（范围为2~sqrt(100000)）后，把这些数的倍数全部筛掉，这样就省去了很多的判断步骤，这种方法称为筛法。这其实就是试商法中的算法原理：若2~sqrt(x)中有整数能够整除x，那么x一定不是素数，也就是若x不是素数，则一定存在某个数位于2~sqrt(x)中（记为i），x是i的倍数。

那么如何实现这个算法呢？我们可以考虑用一个数组a来存储1~N中的所有数，并先把0与1筛去，即初始化数组a，使a[0]=0, a[1]=0, a[2]=2, a[3]=3, ......, a[N]=N；接着对i=2, 3, ......, sqrt(N)分别做：“筛掉a中的所有a[i]的倍数”，即当i+1<=j<=N时，若a[j]是a[i]的倍数，使a[j]=0；最后输出数组中余下的数（a[i]!=0的数）。

据此，我们可编写代码如下：

int i, j, a[N+1] = {0};//假设N为常量
for (i=2; i<=N; i++)
{
    a[i] = i;//初始化a[2]~a[N]
}
for (i=2; i<=sqrt(N); i++)
{
    if (a[i]!=0)//跳过已被筛掉的数
    {
        for (j=i+1; j<=N; j++)
        {
            if (a[j]!=0 && a[j]%a[i]==0)//a[j]!=0的目的是跳过已经被筛掉的数
            {
                a[j] = 0;
            }
        }
    }
}
for (i=0; i<N; i++)
{
    if (a[i])//如果a[i]!=0，则a[i]是素数
    {
        printf("%d\n", a[i]);
    }
}

而在这种方法中，我们发现，当输出素数时，既可以输出i，也可以输出a[i]，那么就说明了有一个值是多余的。其实，我们完全可以使用a[i]当作标记变量，在一开始时把数组a初始化为0，让j遍历2~N，当我们找到非素数则使a[j]变为1，代表筛掉了非素数，最后只需输出元素a[j]=0的下标j即可。据此我们可以优化以上代码：

int i, j, a[N+1] = {0};//假设N为常量
a[0] = 1;
a[1] = 1;//0和1不是素数
for (i=2; i<=sqrt(N); i++)
{
    if (a[i]!=1)//跳过已被筛掉的数
    {
        for (j=i+1; j<=N; j++)
        {
            if (a[j]!=1 && j%i==0)//a[j]!=1的目的是跳过已经被筛掉的数
            {
                a[j] = 1;
            }
        }
    }
}
for (i=0; i<N; i++)
{
    if (!a[i])//如果a[i]==0，则a[i]是素数
    {
        printf("%d\n", i);
    }
}

这样我们就可以省去把数组a初始化为a[2]=2,......,a[N]=N的一步了，从而减少了一个循环，提高了程序运行的效率。

在这里，我们判断j是否为i的倍数时，采用的时让j从i+1遍历到N，判断j对i求余是否为0的方法。这样做其实也是不必要的，因为相比于i的倍数来说，不是i的倍数的数要多得多。因此我们可以让j从2开始，直到i*j>N，这样我们就能直接找到i的倍数i*j，并让a[i*j]=1。

int i, j, a[N+1] = {0};//假设N为常量
a[0] = 1;
a[1] = 1;//0和1不是素数
for (i=2; i<=sqrt(N); i++)
{
    if (a[i]!=1)//跳过已被筛掉的数
    {
        j = 2;
        while (i*j<=N)
        {
            a[i*j] = 1;
            j++;
        }
    }
}
for (i=0; i<N; i++)
{
    if (!a[i])//如果a[i]==0，则a[i]是素数
    {
        printf("%d\n", i);
    }
}

这就是代码经优化的最终版本啦！

---

筛法求素数的核心是依次筛掉2~sqrt(N)中素数的倍数，直到a中仅剩下素数为止，这种方法巧用数组，利用数组作为标志变量，是一种效率较高的求素数的算法。

---

参考文献：

[1]百度百科.