about 4 years ago

埃拉托斯特尼筛法(Sieve of Eratosthenes)

也许这种筛法是大家接触最频繁的筛法吧。也就是每次把当前的质数的倍数全部标记为合数,很好理解,正确性也非常显然。

这种筛法扩展空间非常大,可以拿来求很多东西,但是美中不足的是,每个合数可能会被不止一个质数筛去,因此它的复杂度应该大于\(O(n)\),实际上它的复杂度是\(O(nloglogn)\)。(复杂度证明实在很鬼畜,我也不会,大家也没必要去钻牛角尖研究它,有兴趣的就看看http://www.wikiwand.com/en/Sieve_of_Eratosthenes#/Algorithm_complexity的证明吧)

别小看了那个loglogn,有些题目喜欢把数据范围设得很大,可能就会卡掉这种埃氏筛法。

模板

int primes[MAXN],tot=0;
bool isPrime[MAXN];

void getPrime()
{
    memset(isPrime,true,sizeof(isPrime));
    for(int i=2;i<MAXN;i++)
    {
        if(isPrime[i])
        {
            primes[++tot]=i;
            for(int j=i+i;j<MAXN;j+=i)
                isPrime[j]=false;
        }
    }
}

欧拉筛法(Sieve of Euler)

欧拉筛法的核心思想和埃氏筛法差不多,不过它的复杂度是真正的线性\(O(n)\)!因为它的优化使得每个合数都保证只会被一个质数筛去。

首先引出欧拉筛法筛素数的代码:

模板

int primes[MAXN],tot=0;
bool isPrime[MAXN];

void getPrime()
{
    memset(isPrime,true,sizeof(isPrime));
    for(int i=2;i<MAXN;i++)
    {
        if(isPrime[i])
            primes[++tot]=i;
        for(int j=1;j<=tot;j++)
        {
            if(i*primes[j]>=MAXN) break;
            isPrime[i*primes[j]]=false;
            if(i%primes[j]==0) break;
        }
    }
}

可以发现,欧拉筛法其实和上面的埃氏筛法差不多,但是它枚举出来的i是任意数字,然后用枚举出来的i和已经筛出的质数\(p_j\)组合成一个合数,并筛去这个合数。

但是有一句代码非常特别:if(i%primes[j]==0) break;

正是这一句话,让这个筛法得以拥有线性的时间复杂度。因为它让每个合数都保证只被其最小的质因数筛去。下面看看为什么(感谢http://debug18.com/introduction-to-sieve-method/的讲解):

假设对于合数\(n=p_1^{k_1}...p_m^{k_m}\),这里的i就是外层循环里枚举的i,而\(p_1...p_m\)就是内层枚举的质数

1、所有合数都会被筛去(算法的正确性):\(n\)一定会在i枚举到\(p_1^{k_1-1}p_2^{k_2}...p_m^{k_m}\)时被筛去(此时\(p_1\)一定已经筛出来了,因为它是小于等于i的)

2、所有合数都只会被它最小的质因数筛去(算法的复杂度):若当前有一个质数\(p_x(x>1)\)也筛去了\(n\),那么此时的\(i=p_1^{k_1}...p_x^{k_x-1}...p_m^{k_m}\)。

(1)若\(i>p_x\),那么此时内层循环已经早早break掉了,因为在之前循环里已经枚举到了\(p_1,p_1|i \)

(2)若\(i < p_x\)(显然此时x=m),那这时候\(p_x\)还没被加入进质数表

这样这个筛法的正确性和复杂度都被证明出来了。


线性筛求积性函数(如欧拉函数、莫比乌斯函数、逆元)

之前提到的两种筛法都有很大的拓展空间。下面简单选取欧拉筛法求积性函数的几个经典例子。首先还是讲讲线性筛求一般积性函数\(f(x)\)的过程。

积性函数\(f(x)\)满足:\(f(ab)=f(a)f(b),a,b\)互质。

观察之前欧拉筛法的步骤,筛掉一个合数n的同时,我们可以获得它的最小的质因数p(还记得之前的讲解吗,合数n肯定只会被最小的质因数筛去)。假设p在n中出现了k次,那么我们可以通过\(f(p^k),f(\frac n {p^k})\)得到\(f(n)\)(因为\(p^k和\frac n {p^k} 互质\))。

如果我们在筛去一个合数\(i\)的同时,记录下它里面每个质因数\(p_j\)的次数,那么就能在\(i\)转移到\(i·p_t\)时,在筛去\(i·p_t\)的同时得到它里面任意一个质因数的次数(不过这种做法比较少见,大多数的积性函数可以利用它们的性质,直接从\(f(p^k),f(\frac n {p^k})\)得到\(f(n)\),而不需要记录每个质因数的出现次数)

下面就来介绍一下线性筛求欧拉函数、莫比乌斯函数、逆元的过程。

一、线性筛欧拉函数

欧拉函数\(\phi(n)=n(1-\frac 1 {p_1})(1-\frac 1 {p_2})...(1-\frac 1 {p_k}),p_1...p_k\)是n的k个不同的质因数

那么根据定义式可以得到\(\phi(p)=p-1,p\)是质数

假设当前从\(\phi(i),\phi(p_t)\)转移到\(\phi(ip_t)\),有以下几个式子便于我们筛出欧拉函数:

1、如果\(p_t是在ip_t\)中第一次出现的话(也就是\(p_t\)不整除i),则\(\phi(ip_t)=\phi(i)(p_t-1)\)

2、如果\(p_t不是在ip_t\)中第一次出现的话(也就是\(p_t\)整除i),则\(\phi(ip_t)=\phi(i)p_t\)

式子的推导很简单,大家根据欧拉函数定义式算算\(\phi(i)\)是怎么转移到\(\phi(ip_t)\)就知道了

1和2两种公式的条件判断正好可以借助于线性筛的if(i%primes[j]==0) break;一句话解决,是不是非常巧妙?

下面给出线性筛欧拉函数的代码:

```cpp
int phi[MAXN];
int primes[1100000],tot=0;
bool isPrime[MAXN];

void getPhi()
{
memset(isPrime,true,sizeof(isPrime));
phi[1]=1;
for(int i=2;i {
if(isPrime[i])
{
primes[++tot]=i;
phi[i]=i-1;
}
for(int j=1;j {
if(i*primes[j]>=MAXN) break;
isPrime[i*primes[j]]=false;
if(i%primes[j]==0)
{
phi[i*primes[j]]=phi[i]primes[j];
break;
}
phi[i
primes[j]]=phi[i]*(primes[j]-1);
}
}
}
```

2、线性筛莫比乌斯函数

莫比乌斯函数\(\mu(d)\)的定义如下

(1)若d=1,那么\(\mu(d)=1\)

(2)若d=p_1p_2...p_k(p_1...p_k均为互异质数),那么\(\mu(d)=(-1)^k\)

(3)其他情况下,\(\mu(d)=0\)

莫比乌斯函数最常见的用途便是莫比乌斯反演,除此之外,它还有几个非常好的性质,由于它不在本讲义的讨论范围内,因此具体请参考http://blog.csdn.net/acdreamers/article/details/8542292

假设当前从\(\mu(i),\mu(p_t)\)转移到\(\mu(ip_t)\),有以下几个式子便于我们筛出欧拉函数:

1、如果\(p_t\)是在\(ip_t\)中第一次出现的话(也就是\(p_t\)不整除i),则\(\mu(ip_t)=-\mu(i)\)

2、如果\(p_t\)不是在\(ip_t\)中第一次出现的话(也就是\(p_t\)不整除i),则\(\mu(ip_t)=0\)

这两个转移都非常显然,看看莫比乌斯函数定义式就知道了。

同欧拉函数一样,1和2两种公式的条件判断也可以借助于线性筛的if(i%primes[j]==0) break;一句话解决。

线性筛莫比乌斯函数代码:

```cpp
int miu[MAXN],primes[MAXN],tot=0;
bool isPrime[MAXN];

void getMiu()
{
memset(isPrime,true,sizeof(isPrime));
miu[1]=1;
for(int i=2;i {
if(isPrime[i])
{
miu[i]=-1;
primes[++tot]=i;
}
for(int j=1;j {
if(i*primes[j]>=MAXN) break;
isPrime[i*primes[j]]=false;
if(i%primes[j]==0)
{
miu[i*primes[j]]=0;
break;
}
miu[i*primes[j]]=-miu[i];
}
}
}
```

3、线性求逆元

其实逆元的线性求法和之前提到的两种筛法没什么关系,我更喜欢称这个做法叫逆元的递推求法。但是考虑到贾志鹏的线性筛PPT里提到了逆元的线性求法,我就在这里也说下吧。

逆元其实不是函数,但是我们可以把它看成函数\(f(x)\),逆元(函数)的积性性质也是非常显然的,我就不赘述了。

实际上x和x-p在模p意义下的逆元是相同的,因此我们只需要筛出\(1...p-1\)的逆元就够了。

我们假设模数\(p=ki+r,r < i,1 < i p\)

$$ki+r \equiv 0 \pmod p$$

两边同时乘上\(i^{-1}r^{-1}\):

$$kr^{-1}+i^{-1} \equiv 0 \pmod p$$

$$i^{-1} \equiv -kr^{-1} \pmod p$$

$$i^{-1} \equiv -\lfloor \frac p i\rfloor (p \bmod i)^{-1} \pmod p$$

于是我们可以推出\(i\)的逆元\(f(i)=[-\lfloor \frac p i\rfloor f(p \bmod i)]\bmod p\)

所以逆元可以直接通过递推的方式求出来。

逆元的线性求法代码:

cpp
for(LL i=2;i<MOD;i++)
rev[i]=(MOD-MOD/i)*rev[MOD%i]%MOD;



莫比乌斯反演

相信大家一定对"莫比乌斯反演"这个词耳熟能详。下面我就简单介绍一下莫比乌斯反演公式及它的一个简单应用。

\(F(n),f(n)\)是定义在非负整数集合上的两个函数,并且满足条件\(F(n)=\sum_{d|n}f(d)\),那么我们得到结论:

$$f(n)=\sum_{d|n} \mu (d) F(\frac n d)...①$$

$$f(n)=\sum_{n|d} \mu (\frac d n) F(d)...②$$

有些题目会要求求出f(n)的值,但是往往F(n)比f(n)更好求。如果我们难以求出f(n)的值,但是能找出一个F(n)函数,并且可以很轻松地求出F(n)的值的话,就应该使用莫比乌斯反演。

例题:[BZOJ 2301][HAOI2011]Problem b

题目要求的是\(\sum_{a\leq x \leq b}\sum_{c\leq y \leq d}[gcd(x,y)=k]\)的值。

实际上我们可以把问题转化为求\(\sum_{1\leq x \leq n}\sum_{1\leq y \leq m}[gcd(x,y)=k]\),然后做四次计算,用容斥原理求出最终答案。

考虑已知\(n,m\),上面的式子怎么求出来。

我们不妨设一个函数f(k):

$$f(n,m,k)=\sum_{1\leq x \leq n}\sum_{1\leq y \leq m}[gcd(x,y)=k]$$

$$=\sum_{1\leq x \leq \lfloor \frac n k \rfloor}\sum_{1\leq y \leq \lfloor \frac m k \rfloor}[gcd(x,y)=1]$$

$$=f(\lfloor \frac n k \rfloor,\lfloor \frac m k \rfloor,1)$$

f(k)就是x在[1,n],y在[1,m]中的gcd(x,y)=k的有序数对(x,y)个数。

我们再设一个函数F(k):

$$F(n,m,k)=\sum_{1\leq x \leq n}\sum_{1\leq y \leq m}[gcd(x,y)\bmod k=0]=\lfloor \frac n k \rfloor\lfloor \frac m k \rfloor$$

F(k)就是x在[1,n],y在[1,m]中的gcd(x,y)为k的倍数的有序数对(x,y)个数。

显然F(k)和f(k)两个函数满足莫比乌斯反演的条件。则:

$$f(n,m,k)=f(\lfloor \frac n k \rfloor,\lfloor \frac m k \rfloor,1)$$

$$=\sum_{1|n,n\leq \min(\lfloor \frac n k \rfloor,\lfloor \frac m k \rfloor)} \mu (\frac n 1) F(\lfloor \frac n k \rfloor,\lfloor \frac m k \rfloor,n)$$

$$=\sum_{1\leq n\leq \min(\lfloor \frac n k \rfloor,\lfloor \frac m k \rfloor)} \mu (n) F(\lfloor \frac n k \rfloor,\lfloor \frac m k \rfloor,n)$$

我们利用分块+前缀和优化可以将时间复杂度降为\(O(\sqrt {\max(n,m)})\)

补充

贾志鹏的线性筛PPT非常不错,值得一读

← [NOIP 2015复习]图论篇 [NOIP 2015复习]图论篇 →
 
comments powered by Disqus