数论——容斥原理、莫比乌斯函数

1、容斥原理：

• 时间复杂度为$O(2^N)$，下面会有证明。
• 举一个简单的例子：用韦恩图来思考，求$S1$、$S2$、$S3$三个集合的原有元素的并集，那么结果为：$S1+S2+S3-S1 \cap S2-S1 \cap S3-S2 \cap S3+S1 \cap S2 \cap S3$。
  
• 以此类推到$N$个圆的交集：用容斥原理的方法答案为所有单个集合的元素个数-所有两个集合互相交集的元素个数+所有三个集合互相交集的元素个数…
• 我们知道容斥原理公式一共涉及到的元素个数为：$C_N^1+C_N^2+C_N^3+…+C_N^N$。因为$C_N^0+C_N^1+C_N^2+C_N^3+…+C_N^N=2^n$，因此$C_N^1+C_N^2+C_N^3+…+C_N^N=2^n-1$，因此容斥原理公式一共涉及到的元素个数为$2^n-1$。关于此公式($C_N^0+C_N^1+C_N^2+C_N^3+…+C_N^N=2^n$)的证明，我们可以假设等号左边为对于$N$个物品所有选法的总个数，等号右边考虑每个物品选与不选两种情况，因此等式成立。
• 因此容斥原理的时间复杂度为$O(2^N)$。
• 容斥原理的证明：对于容斥原理$|S1 \cup S2 \cup … \cup SN|=\sum_{i=1}^N{Si}-\sum_{i, j}^N{Si \cap Sj}+\sum_{i, j, k}^N{Si \cap Sj \cap Sk}+…$
  对于一个元素$x$，它在$k$个集合中，$1 \leq k \leq N$，它本身被选择的次数为$C_k^1-C_k^2+C_k^3-…+(-1)^{k-1}C_k^k$。我们知道一个结论：$C_k^1-C_k^2+C_k^3-…+(-1)^{k-1}C_k^k=1$，因此对于每一个元素$x$，它只被计算了$1$次，证毕。

例题：AcWing 890. 能被整除的数

给定一个整数$n$和$m$个不同的质数$p1,p2,…,pm$。请你求出$1$到$n$中能被$p1,p2,…,pm$中的至少一个数整除的整数有多少个。

首先我们知道，在$N$个数中能被$x$整除的数的个数为$\lfloor{N/x}\rfloor$。

因此我们只需要根据容斥原理，求出可以被单个元素整除的个数之和-可以被两个元素整除的个数之和+可以被三个元素整除的个数之和…我们用位运算来求得答案，时间复杂度为$O(2^N)$。

#include <bits/stdc++.h>
#define int long long
using namespace std;
int p[20];

void work() {
    int n, m;
    cin >> n >> m;
    for (int i = 0; i < m; i++) cin >> p[i];

    int res = 0;
    for (int i = 1; i < 1 << m; i++) {
        int t = 1, s = 0;
        for (int j = 0; j < m; j++)
            if (i >> j & 1) {
                if (t * p[j] > n) {
                    t = -1;
                    break;
                } 
                t *= p[j];
                s++;
            }
        if (t != -1) {
            if (s % 2) res += n / t;
            else res -= n / t;
        }
    }
    cout << res << endl;
}

int32_t main() {
    ios::sync_with_stdio(false);
    cin.tie(0);
    cout.tie(0);
    int cas;
    cas = 1;
    while (cas--) work();
}

例题：AcWing 214. Devu和鲜花

有$N$个盒子，第$i$个盒子中有$Ai$枝花。同一个盒子内的花颜色相同，不同盒子内的花颜色不同。要从这些盒子中选出$M$枝花组成一束，求共有多少种方案。若两束花每种颜色的花的数量都相同，则认为这两束花是相同的方案。结果需对$10^9+7$取模之后方可输出。

我们先考虑：从$N$个盒子选$M$枝花，每个盒子的花的个数为无限个，问一共能选多少枝花？

此问题等价于从$N$个盒子选$M+N枝$花，那么每个盒子至少选$1$枝。那么此问题由等价于把$N+M$个点分成$N$份，我们可以用隔板法来做，一共有$N+M-1$个空隙，有$N-1$个板子，因此答案为$C_{N+M-1}^{N-1}$。

拓展到此问题，第$i$个盒子中有$Ai$枝花。那么我们可以反过来考虑，用总共的答案$C_{N+M-1}^{N-1}$减去其中第$i$个盒子被拿走了大于$Ai$枝花的方案。第$i$个盒子被拿走了大于$Ai$枝花的方案数为：假设此盒子已经被拿走了$Ai+1$枝花，那么等价于前面的问题，从$N$个盒子中共拿走$M-Ai-1$枝花的方案数，等价于从N个盒子拿走$M-Ai-1+N$的方案数，每个盒子至少被拿$1$枝，因此答案为$C_{M-Ai-1+N-1}^{N-1}$。

根据容斥原理来做，可知答案为总共的$C_{N+M-1}^{N-1}$减去所有$1$个盒子不满足的加上所有$2$个盒子不满足的减去所有$3$个盒子不满足的…

#include <bits/stdc++.h>
#define int long long
using namespace std;
int A[20];
constexpr int mod = 1e9 + 7;
int down = 1;

int qmi(int a, int b, int p) {
    int res = 1;
    while (b) {
        if (b & 1) res = res * a % p;
        a = a * a % p;
        b >>= 1;
    }
    return res;
}

int C(int a, int b) {
    if (a < b) return 0;
    int up = 1;
    for (int i = a; i > a - b; i--) up = i % mod * up % mod;
    return up * down % mod;
}

void work() {
    int n, m;
    cin >> n >> m;
    for (int i = 0; i < n; i++) cin >> A[i];
    for (int j = 1; j <= n - 1; j++) down = j * down % mod;
    down = qmi(down, mod - 2, mod);

    int res = 0;
    for (int i = 0; i < 1 << n; i++) {
        int a = n + m - 1, b = n - 1;
        int sign = 1;
        for (int j = 0; j < n; j++)
            if (i >> j & 1) {
                sign *= -1;
                a -= A[j] + 1;
            }
        res = (res + C(a, b) * sign) % mod;
    }
    cout << (res % mod + mod) % mod << endl;
}

int32_t main() {
    ios::sync_with_stdio(false);
    cin.tie(0);
    cout.tie(0);
    int cas;
    cas = 1;
    while (cas--) work();
}

2、莫比乌斯函数：

我们举例一道经典的应用题，求$1$到$N$中与$a$互质的数的个数：那么根据容斥原理，设$S_i$为$1$到$N$中和$a$有公因子i的数的个数，答案为$N-S_2-S_3-S_5-S_7..+S_{2,3}+S_{3,5}+S_{2,5}+…$，我们可以惊奇的发现，其答案为$\sum_{i=0}^{N}{u(i)*S_i}$。

我们可以根据线性筛质数在$O(N)$的时间内算出前$N$个数的莫比乌斯数。

int primes[N], cnt;
bool st[N];
int mobius[N];

void init(int n) {
    mobius[1] = 1;
    for (int i = 2; i <= n; i++) {
        if (!st[i]) {
            primes[cnt++] = i;
            mobius[i] = -1;
        }
        for (int j = 0; primes[j] * i <= n; j++) {
            int t = primes[j] * i;
            st[t] = true;
            if (i % primes[j] == 0) {
                mobius[t] = 0;
                break;
            }
            mobius[t] = mobius[i] * -1;
        }
    }
}

AcWing 215. 破译密码

对于给定的整数$a$，$b$和$d$，有多少正整数对$x$，$y$，满足$x \leq a$，$y\leq b$，并且$gcd(x, y)=d$。$5e4$组询问，$a、b$的数据范围为$5e4$。

根据数据的范围，我们可以推断出时间复杂度为$O(nlogn)$。

每次询问的问题等价于：有多少正整数对$x$，$y$，满足$x \leq a/d$，$y\leq b/d$，并且$gcd(x, y)=1$。那么根据容斥原理：值为$min(a/d, b/d)-S_2-S_3-S_5+S_{2,3}+S_{2,5}+S_{3,5}-S{2,3,5}…$，其答案为$\sum_{i=0}^{min(a/d, b/d)}{u(i)*S_i}$，也就是$\sum^{min(a,b)}_i=a/i∗b/i∗u[i]$。我们可以推断出时间复杂度为$O(n)$。

考虑把上述过程优化，发现，这个式子中虽然i要枚举$N$次，但是实际上因为整除的原因$ai$的值很少，只有$2\sqrt a$个。

因为$a/1 、a/2、a/3、…$是单调递减的，并且有的值相同，所以整个序列一共有有$2\sqrt a$个值。（证明：在分母为$1$到$\sqrt a$之间，值的个数为$a / \sqrt a$个值，在$\sqrt a +1$到$n$之间，值的个数为$a / \sqrt a$个值）。

设$g(x)$表示$a/x$的取值不变的最大的$x$值，那么$a/x=a/g(x)$，并且$a/x>a/(g(x)+1)$，其中$g(x)=a/(a/x)$。

证明$a/x=a/g(x)$：

证明：$g(x)=a/(a/x)$：

综上：将原来的序列分成$2\sqrt a$段，而且每次都会跳一段，所以总共会跳$2\sqrt a$次，时间复杂度就是$O(\sqrt a)$。

加上询问后总的时间复杂度为$O(nlogn)$。

#include <bits/stdc++.h>
typedef long long ll;
using namespace std;
const int N = 50010;

int primes[N], cnt;
bool st[N];
int mobius[N], sum[N];

void init(int n) {
    mobius[1] = 1;
    for (int i = 2; i <= n; i++) {
        if (!st[i]) {
            primes[cnt++] = i;
            mobius[i] = -1;
        }
        for (int j = 0; primes[j] * i <= n; j++) {
            int t = primes[j] * i;
            st[t] = true;
            if (i % primes[j] == 0) {
                mobius[t] = 0;
                break;
            }
            mobius[t] = mobius[i] * -1;
        }
    }
    for (int i = 1; i <= n; i++) sum[i] = sum[i - 1] + mobius[i];
}

void work() {
    int a, b, d;
    cin >> a >> b >> d;
    a /= d, b /= d;
    int n = min(a, b);
    ll res = 0;
    for (int l = 1, r; l <= n; l = r + 1) {
        r = min(n, min(a / (a / l), b / (b / l)));
        res += (sum[r] - sum[l - 1]) * (ll)(a / l) * (b / l);
    }
    cout << res << endl;
}

int32_t main() {
    ios::sync_with_stdio(false);
    cin.tie(0);
    cout.tie(0);
    init(N - 1);
    int cas;
    cin >> cas;
    while (cas--) work();
}