线性筛与莫比乌斯反演

目录

线性筛与莫比乌斯反演

和上篇文章一样，一直没有研究这个东西，结果又考了GG……TAT
下定决心学一学，搞好这个东西。

线性筛

筛质数有很多方法，好像很厉害的有洲阁筛 \(O(\frac{n^{\frac{3}{4}}}{\log n})\) 、杜教筛 \(O(n^{\frac{2}{3}})\)（然而我都不会QAQ）(其实这些不是用来筛素数的2333 用来筛积性函数前缀和的 \(Update \ 2018 \cdot 4 \cdot 9\) )，还有暴力筛（就是枚举一个数的倍数）复杂度是 \(O(n \ln n)\) 的。
我只学了比较简单而且实用的线性筛法。
这种筛法是避免一个数被重复筛几遍，所以效率均摊下来可以达到线性。(网上有证明)

代码

const int N = 100000;
bool is_prime[N+100];
int prime[N], cnt = 0;
void find_prime() {
    Set(is_prime, true);
    is_prime[0] = is_prime[1] = false;
    For(i, 2, N) {
        if (is_prime[i]) prime[++cnt] = i;
	For(j, 1, cnt) {
		if (i * prime[j] > N) break;
		is_prime[i * prime[j]] = false;
		if (i % prime[j] == 0) break; //here
	}
    }
}

讲解

这个代码有一个关键点 就是上面的\(here\) 这个意义就是对于一个合数\(m\)可以分解为\(m=p_1^{r_1}*...*p_n^{r_n}\)其中
\(p_i\)为质数，那么我们筛\(m\)的时候之前把\(p_1\)筛掉了，所以在枚举\(i\)的时候。

  1. 如果\(i\)为素数没问题，直接向后继续推（因为筛出的质数都类似\(m=p_1*p_2\)的形式，所以不可能重复）。
  2. 如果为合数，那么\(i\)可以分解成\(i={p_1}^{r_1}*...*{p_n}^{r_n}\)形式其中\(p_1-p_n\)是递增的，
    那么\(p_1\)是最小的那个质数。\(i \bmod p_1 = 0\)的时候，就不用继续枚举了，所以我们就只能筛出不大于\(p_1\)的质数\(*i\)。

复杂度好像是线性的（ \(O(n)\) ），我不太会证复杂度。。

莫比乌斯反演

莫比乌斯反演很多时候都能大大简化运算……

定理

\(F(n)\)和\(f(n)\)是定义在非负整数集合上的两个函数，并且满足条件$$F(n)=\sum \limits \limits _{d|n}{f(d)}$$。那么我们就能得到结论：

\[f(n)=\sum_{d|n}\mu(d)F(\frac{n}{d}) \]

在上面的公式中有一个\(\mu(d)\)函数（莫比乌斯函数），它的定义如下：

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

  2. 若\(d=p_1p_2...p_k\)，\(p_i\)均为互异质数，那么\(\mu(d)=(-1)^{k}\)。这个我的理解就是\(d\)的质因数个数为偶数的话，那么\(\mu(d)=1\) 否则为 \(-1\) 。

  3. 其他情况下\(\mu(d)=0\)这个就是对上面那条的拓展了，就是指的\(d\)没有一个平方因子，或者说没有一个质因子的次数大于 \(1\) 。

代码

const int N = 100100;
bool is_prime[N+100];
int mu[N+100] = {0, 1}, cnt = 0, prime[N+100];
void init() {
	Set(is_prime, true);
	is_prime[1] = false;
	For (i, 2, N) {
		if (is_prime[i]) {
			prime[++cnt] = i;
			mu[i] = -1; //质数的质因子个数肯定为奇数个就是1
		}
		For (j, 1, cnt) {
			if (i * prime[j] > N) break;
			is_prime[i * prime[j]] = false;
			if (i % prime[j]) mu[i * prime[j]] = -mu[i]; //多了一个质因子直接变为原来结果的相反数
			else {
				mu[i * prime[j]] = 0; //这个将要被筛的数至少具有两个prime[j]的因子
				break;
			}
		}
	}
}

常见的定理

然后还要提一下的就是一些常见的定理，证明嘛……

一般都是先分解质因数，然后再根据组合数性质去算，比如第一个。

要么就是对于一些常见的反演格式进行反演，比如第二个。

\[\sum _{d|n} \mu(d)=[n=1] \]

\[\sum _{d|n} \frac{\mu(d)}{d}=\frac{\varphi(n)}{n} \]

莫比乌斯反演的证明

证明：

\[\begin{aligned} \sum_{d|n}\mu(d)F(\frac{n}{d})&=\sum_{d|n}\mu(d)\sum \limits_{d'|\frac{n}{d}}f(d')\\ &=\sum_{d'|n}f(d')\sum_{d|\frac{n}{d'}}\mu(d)\\&=f(n) \end{aligned} \]

Q.E.D

一些例题（难题）

Luogu 【P1829】[国家集训队]Crash的数字表格

题意

\[\sum_{i=1}^{n} \sum_{i=1}^{m} \mathrm{lcm}(i,j) \ (n,m \le 10^7) \]

题解

一个莫比乌斯反演然后化式子。

\[\begin{aligned} ans &= \sum _{i=1}^{n} \sum _{i=1}^{m} \mathrm{lcm}(i,j)\\ &= \sum _{i=1}^{n} \sum _{i=1}^{m} \frac{ij}{\gcd(i,j)}\\ &= \sum _{d=1}^{min(n,m)} d \sum _{i=1}^{\lfloor \frac{n}{d} \rfloor} \sum _{j=1}^{\lfloor \frac{m}{d} \rfloor} ij \ [\gcd(i,j)=1]\\ \end{aligned} \]

这个就是一个更换枚举相的操作了，是个套路。你先枚举所有可能的\(\gcd\)再计算这种\(\gcd\)的贡献。

比如前面的那个\(d\)就是我们枚举的\(\gcd\)，后面所有可能的数对，就是在\(\lfloor \frac{n}{d} \rfloor\)和\(\lfloor \frac{m}{d} \rfloor\)中的所有互质的数对的乘积在乘上\(d\)。

这个可以简单理解一下，就是两个数分别除以他们的最大公因数，然后两个数肯定是互质的。
但其对于答案的贡献就多除以了一个\(d\)，所以要乘回来。

\[\begin{aligned} ans =\sum_ {d=1}^{\min(n,m)} d \sum_{i=1}^{\lfloor \frac{n}{d} \rfloor} \sum_{j=1}^{\lfloor \frac{m}{d} \rfloor} \sum_{x|\gcd(i,j)}\mu(x) * i * j \end{aligned} \]

这个就是运用了前面的公式\(\sum \limits _{d|n} \mu(d)=[n=1]\)来替代了\([gcd(i,j)=1]\)的条件。（这个就是套路了）

然后我们继续推：

\[\begin{aligned} ans &=\sum_{d=1}^{\min(n,m)} d \sum_{x=1}^{min(\lfloor \frac{n}{d} \rfloor,\lfloor \frac{m}{d} \rfloor)}\mu(x) \ x^2 \sum_{i=1}^{\lfloor \frac{n}{dx} \rfloor} i \sum_{j=1}^{\lfloor \frac{m}{dx} \rfloor} j \end{aligned} \]

这个也是套路，把\(x\)提前了。就是改成了枚举\(x\)看看它的对于答案的贡献是多少。

很容易发现，就是在\([1,\lfloor \frac{n}{dx} \rfloor]\)中的所有数乘上\(x\)

就是原来可行的\(i\)。然后我们就可以根据这个来优化了。

前面那两个\(\sum \limits\)就是\(O(n \ln \ n)\)（令\(n=\min(n,m)\)）的复杂度。（就是调和级数 \(H(n) = \sum_{i=1}^{n} \frac{n}{i}\)）。后面的那两个，直接用等差数列求和公式\(O(1)\)算。

但这个仍然过不去...（\(O(1.66*10^8)\)，\(\bmod\)的常数还很大）所以就需要来用套路的整除分块了。
就是把后面两个 \(\sum \limits\) 很多一样答案的地方一起处理掉，所以对于那个 \(\mu(x) \ x^2\) 还要记一个前缀和。

总复杂度\(O(\sum \limits _{i=1}^{n} \sqrt{\frac{n}{i}})=O(pass)\)这个我也不会算。。会算的大佬私聊啊23333

代码

#include <bits/stdc++.h>
#define For(i, l, r) for (int i = (l), _end_ = (int)(r); i <= _end_; ++i)
#define Fordown(i, r, l) for (int i = (r), _end_ = (int)(l); i >= _end_; --i)
#define Set(a, v) memset(a, v, sizeof(a))
using namespace std;

bool chkmin(int &a, int b) { return b < a ? a = b, 1 : 0; }
bool chkmax(int &a, int b) { return b > a ? a = b, 1 : 0; }

void File() {
#ifdef zjp_shadow
    freopen("P2154.in", "r", stdin);
    freopen("P2154.out", "w", stdout);
#endif
}

typedef long long ll;
ll n, m;
const int N = 1e7 + 1e3;
const ll Mod = 20101009, inv2 = (Mod + 1) / 2;
bool is_prime[N];
int mu[N], prime[N], cnt;
ll sum[N];

inline void add(ll &x, ll y) {
    x = ((x + y) % Mod + Mod) % Mod;
}

void Get_Mu(int maxn) {
    int res;
    Set(is_prime, true);
    is_prime[0] = is_prime[1] = false;
    mu[1] = 1;
    For(i, 2, maxn) {
        if (is_prime[i])
            prime[++cnt] = i, mu[i] = -1;
        For(j, 1, cnt) {
            if ((res = i * prime[j]) > maxn) break;
            is_prime[res] = false;
            if (i % prime[j]) mu[res] = -mu[i];
            else { mu[res] = 0; break; }
        }
    }
    For(i, 1, maxn) add(sum[i], sum[i - 1] + (ll)mu[i] * i * i % Mod);
}

inline ll fsum(ll x) { return x * (1 + x) % Mod * inv2 % Mod; }

int main() {
    File();
    cin >> n >> m;
    Get_Mu(min(n, m));
    ll ans = 0;
    For(d, 1, min(n, m)) {
        int n_ = n / d, m_ = m / d;
        For(x, 1, min(n_, m_)) {
            if (!mu[x]) continue;
            int Next = min(n_ / (n_ / x), m_ / (m_ / x));
            add(ans, 1ll * d * (sum[Next] - sum[x - 1]) % Mod * fsum(n_ / x) % Mod * fsum(m_ / x) % Mod);
            x = Next;
        }
    }
    cout << ans << endl;
    return 0;
}

BZOJ3994 [SDOI2015]约数个数和

题意

求$$\sum \limits _{i=1}^{n} \sum \limits _{j=1}^{m} d(ij)$$，\(d(x)\)表示\(x\)的约数个数。(\(n,m \le 10^5\))

题解

一个反演。当初数学一本通没看懂，真是本垃圾书。

肖大佬博客一看就懂了2333。只是中间有一步化\([gcd(i,j)=1]\)还是习惯变两步，容易理解些QwQ

posted @ 2018-01-23 21:45  zjp_shadow  阅读(1681)  评论(3)  编辑  收藏  举报