数论基础

此笔记是我学习 ACwing算法基础课 的学习笔记。内容很多，请耐心学习，可以去网站上提交测试代码。

数论

质数（素数）

在大于1的整数中，如果只包含1和本身这两个约数，就被称为质数或素数

质数的判定——试除法

题目描述

给定n个正整数$a_i$，判定每个数是否是质数。

输入格式:

第一行包含整数n。

接下来n行，每行包含一个正整数$a_i$。

输出格式

共n行，其中第 i 行输出第 i 个正整数$a_i$是否为质数，是则输出“Yes”，否则输出“No”。

数据范围

$1≤n≤100$,
$1≤a_i≤2∗109$

输入样例：

2
2
6

输出样例：

Yes
No

实现代码:

#include <iostream>

using namespace std;

bool is_prime(int n) {
    if (n < 2) return false;
    for (int i = 2; i < n; ++i) if (n%i == 0) return false;
        return true;
}

int main() {
    int n, t;
    cin >> n;
    for (int i = 0; i < n; ++i) {
        cin >> t;
        if (is_prime(t)) cout << "Yes\n";
        else cout << "No\n";
    }
    return 0;
}

这时候的算法时间复杂度为$O(n)$。

但是有一个性质：$若d|n，则\frac{n}{d}|n$，例如$2|12，则6|12$，因此在判断12是否是质数时，只需要枚举一个数（2），即可，假设$d \leq \frac{n}{d}$,于是循环即可修改为：

for (int i = 2; i <= n/i; ++i) if (n%i == 0) return false;

这时候的算法时间复杂度为$O(\sqrt n)$

分解质因数——试除法

题目描述

给定n个整数$a_i$，将每个数分解质因数，并按照质因数从小到大的顺序输出每个质因数的底数和指数

输入格式

第一行包含整数n。

接下来n行，每行包含一个整数$a_i$

输出格式

对于每个正整数$a_i$，按照从小到大的顺序输出其分解质因数后，每个质因数的底数和指数，每个底数和指数占一行。

每个正整数的质因数全部输出完毕后，输出一个空行。

数据范围

$1 \leq n \leq 100$,

$1 \leq a_i \leq 2*10^9$

输入样例：

2
6
8

输出样例：

2 1
3 1

2 3

说明：$6 = 2^1 * 3^1;8 = 2^3$

实现代码:

#include <iostream>

using namespace std;

void divide(int n) {
    for (int i = 2; i <= n; ++i)
        if (n%i == 0) {
            int s = 0;
            while (n%i == 0) {
                n /= i;
                s++;
            }
            cout << i << " " << s << endl;
        }

    if (n > 1) cout << n << " 1\n";
    cout << endl;
}

int main() {
    int n, t;
    cin >> n;
    for (int i = 0; i < n; ++i) {
        cin >> t;
        divide(t);
    }
    return 0;
}

该方法依次判断每一个数是否能被模尽（n%i == 0），如果能，则就是答案

这里看起来好像是每一个数都在判断，看似包含了合数。但实际上我们知道，若$i|n \and j|i$则$j|n \and j<i$。因此遇到的第一个可以模尽的数，它应该是没有约数的。即当n%i == 0成立时，i没有其他约数，因此i就是质数。

优化：

根据性质：n中最多只包含一个大于$\sqrt{n}$的质因子。于是可以把循环条件改为：

for (int i = 2; i <= n/i; ++i)

因此时间复杂度为：$O(\log n)到O(\sqrt{n})$

筛质数

其思想就是从2开始筛选（1一定不是质数），对于一个数，把它的所有在范围内的倍数筛选掉（因为它的倍数一定是非质数）。当遍历到一个数x时，如果它没有被筛选掉，则说明在2～x-1中没有它的约数（在x-1～N中更不可能有它的约数），因此它是质数，则添加到质数数组中。

题目描述

给定一个正整数n，请你求出1～n中质数的个数

输入格式

共一行，包含整数n

数据范围

$1 \le n \le 10^6$

输入样例

8

输出样例

4

实现代码

埃式筛法：

#include <iostream>

using namespace std;

const int N = 1e6+1;
int primes[N], cnt; // primes存质因数，cnt表示个数
bool st[N]; // 表示是否筛掉

void get_prime(int n) {
    for (int i = 2; i <= n; ++i) {
        if (!st[i]) primes[cnt++] = i;
        for (int j = i+1; j <= n; j+=i) st[j] = true;
    }
}

int main() {
    int n, t;
    cin >> n;
    get_prime(n);
    cout << cnt << endl;
    return 0;
}

这种筛法就是直接根据上诉思想进行筛选，它的时间复杂度为$O(n\log n)$

线性筛法：

void get_prime(int n) {
    for (int i = 2; i <= n; ++i) {
        if (!st[i]) primes[cnt++] = i;
        for (int j = 0; primes[j] <= n/i; ++j) {
            st[primes[j] * i] = true;
            if (i%primes[j] == 0) break; // primes[j] 是i的最小质因子
        }
    }
}

这种筛法特点是只会被最小质因子筛掉，并且每个合数只会被筛掉一次。

原因是它不是无脑筛选i的每一个倍数（埃式筛对每个数i都会筛选其$1、2、3、\cdots、k$），而是筛选掉primes[j]*i，即已有质数倍，i不同则primes[j]*i一定不同。

if (i%primes[j] == 0) break;当遇到了i的第一个质因数的时候，就停止筛选，避免重复筛选$primes[j]^2i$(到后面这个数一定会被$primes[j]i$筛选一次的)，同理，避免重复筛选$primes[j]^3i、primes[j]^4i、\cdots、primes[j]^k*i$。这样可以保证同一个合数不会被筛选两次。

直观一点：1 2 3 4 5 6 7 8 9 10 11 12

• 当i=2时：

  primes[] = [2]

  筛选掉2*2=4

• 当i=3时：

  primes[] = [2, 3]

  筛选掉2*3=6,3*3=9

• 当i=4时：

  primes[] = [2, 3]

  筛选掉：2*4=8

  然后4%2==0停止筛选

• 当i=5时：

  primes[] = [2, 3, 5]

  筛选掉：2*5=10(3 > 12/5，因此3、5与5的倍数不用再筛选)

• 当i=6时：

  primes[] = [2, 3, 5]

  筛选掉：2*6=12

  举例就到这里了，可以发现每一个数都只被筛选了一次。其中被筛选掉的12是比较特殊的，它的质因数为2、3，且有算数基本定理可以表达为$2^2*3^1$，因此12可以被4筛选一次，也可以被6筛选一次，所以为了避免麻烦，就让6筛选掉即可。

这个方法更快。

约数

试除法求约数

题目描述

给定n个正整数$a_i$,对于每个整数$a_i$，请你按照从小到大的顺序输出它的所有约数。

输入格式

第一行包含整数n。

接下来n行，每行包含一个整数$a_i$。

输出格式

输出共n行，其中第i行输出第i个整数$a_i$的所有约数。

数据范围

$1 \le n \le 100$,

$2 \le a_i \le 2*10^9$

输入样例

2
6
8

输出样例

1 2 3 6
1 2 4 8

实现代码

#include <iostream>
#include <vector>

using namespace std;

vector<int> get_divisors(int n) {
    vector<int> res;
    for (int i = 1; i <= n/i; ++i) {
        if (n%i == 0) res.push_back(i);
        if (i != n/i) res.push_back(n/i);
    }
    sort(res.begin(), res.end());
    return res;
}

int main() {
    int n, t;
    cin >> n;
    for (int i = 0; i < n; ++i) {
        cin >> t;
        auto res = get_divisors(t);
        for (auto x: res) cout << x << " ";
        cout << endl;
    }
    return 0;
}

时间复杂度为$O(\sqrt n)$

约数个数

性质

由算数基本定理：$N=p_1^{\alpha_1} + p_2^{\alpha_2} + \cdots + p_k^{\alpha_k}$，可以得出，约数个数为$(\alpha_1+1)(\alpha_2+1)\cdots(\alpha_k+1)$。

算术基本定理，又称为正整数的唯一分解定理，即：每个大于1的自然数，要么本身就是质数，要么可以写为2个或以上的质数的积，而且这些质因子按大小排列之后，写法仅有一种方式。

例如：$6936=2^3×3×17^2，1200=2^4×3×5^2，5207=41×127$。

题目描述

给定n个正整数$a_i$,请你输出这些数的乘积的约数个数，答案对$10^9+7$取模

输入格式

第一行包含整数n。

接下来n行，每行包含一个整数$a_i$。

输出格式

输出一个整数，表示所给正整数的乘积的约数个数，答案需对$10^9+7$取模。

数据范围

$1 \le n \le 100$,

$2 \le a_i \le 2*10^9$

输入样例

3
2
6
8

输出样例

12

实现代码

#include <iostream>
#include <unordered_map>

using namespace std;

typedef long long ll;
const int mod = 1e9+7;
unordered_map<int, int> primes;
ll res = 1;

void divide(int n) {
    for (int i = 2; i <= n; ++i)
        while (n%i == 0) {
            n /= i;
            primes[i] ++;
        }

    if (n > 1) primes[n] ++;
}

int main() {
    int n, t;
    cin >> n;
    for (int i = 0; i < n; ++i) {
        cin >> t;
        divide(t);
    }
    
    for (auto prime: primes)
        res = res*(prime.second+1)%mod;
    cout << res << endl;
    return 0;
}

这个代码利用到了分解质因数——试除法的内容，来求底数$p_i$和指数$\alpha_i$（利用哈希表来存储）。然后利用性质，进行计算。

约数之和

性质

由算数基本定理：$N=p_1^{\alpha_1} + p_2^{\alpha_2} + \cdots + p_k^{\alpha_k}$，可以得出，约数之和为$(p_1^0+p_1^1+\cdots+p_1^{\alpha_1})\cdots(p_k^0+p_k^1+\cdots+p_k^{\alpha_k})$

题目描述

给定n个正整数$a_i$,请你输出这些数的乘积的约数之和，答案对$10^9+7$取模

输入格式

第一行包含整数n。

接下来n行，每行包含一个整数$a_i$。

输出格式

输出一个整数，表示所给正整数的乘积的约数之和，答案需对$10^9+7$取模。

数据范围

$1 \le n \le 100$,

$2 \le a_i \le 2*10^9$

输入样例

3
2
6
8

输出样例

252

实现代码

#include <iostream>
#include <unordered_map>

using namespace std;

typedef long long ll;
const int mod = 1e9+7;
unordered_map<int, int> primes;
ll res = 1;

void divide(int n) {
    for (int i = 2; i <= n; ++i)
        while (n%i == 0) {
            n /= i;
            primes[i] ++;
        }

    if (n > 1) primes[n] ++;
}

int main() {
    int n, t;
    cin >> n;
    for (int i = 0; i < n; ++i) {
        cin >> t;
        divide(t);
    }
    
    for (auto prime: primes) {
        int p = prime.first, a = prime.second;
        ll tp = 1;
        while (a--) tp = (tp*p+1)%mod;
        res = res * tp %mod;
    }
    cout << res << endl;
    return 0;
}

这个代码利用到了分解质因数——试除法的内容，来求底数$p_i$和指数$\alpha_i$（利用哈希表来存储）。然后利用性质，进行计算。

最大公约数(欧几里得算法、辗转相除法)

题目描述

给定n对正整数$a_i,b_i$,请你求出每对数的最大公约数。

输入格式

第一行包含整数n。

接下来n行，每行包含一个整数对$a_i,b_i$。

输出格式

输出共n行，每行输出一个整数对的最大公约数。

数据范围

$1 \le n \le 10^5$,

$1 \le a_i,b_i \le 2*10^9$

输入样例

2
3 6
4 6

输出样例

3
2

实现代码

#include <iostream>

using namespace std;

int gcd(int a, int b) {
    return b ? gcd(b, a%b) : a;
}

int main() {
    int n, a, b;
    cin >> n;
    while (n--) {
        cin >> a >> b;
        cout << gcd(a, b) << endl;
    }
}

时间复杂度为$O(\log n)$

欧拉函数

1～N中与N互质的个数被称为欧拉函数，记为$\varPhi(N)$。

若在算数基本定理中，$N=p_1^{\alpha_1}p_2^{\alpha_2}\cdots p_m^{\alpha_m}$,则：

$\varPhi(N) = N\frac{p_1-1}{p_1}\frac{p_2-1}{p_2}\cdots\frac{p_m-1}{p_m} = N(1-\frac{1}{p_1})(1-\frac{1}{p_2})\cdots (1-\frac{1}{p_m})，\varPhi(1) = 1$

例如：

$1,2,3,4,5,6$,其中，只有1，5与6互质，因此：

$\varPhi(6) = 2$​​​

用算数基本定理计算：

$6=2^13^1$，$\varPhi(6) = 6\frac{2-1}{2}\frac{3-1}{3} = 6\frac{1}{2}*\frac{2}{3} = 2$

---

互质：若N个整数的最大公因数是1，则称这N个整数互质（$N \ge 2$​）。

推导（根据容斥定理）：

1. 从 1～N 中去掉$p_1,p_2,\cdots,p_k$的所有倍数
2. 加上所有$p_i*p_j$的倍数（加上交集）
3. 减去所有$p_ip_jp_k$的倍数（多加了，又要再减一次交集）
4. 加上所有$p_ip_jp_k$的倍数（上一步减交集后，因为多加了三次交集，又减去了三次交集，因此还需要加上一次），减去所有$p_ip_jp_k*p_l$的倍数（因为又会有交集）
5. 以次类推，就得到了上述式子

容斥定理：

以三个集合举例：

就是要求元素个数，就是$A \cup B \cup C - A\cap B - B\cap C-A \cap C + A \cap B\cap C$​

如果推广到K个集合，就是上诉式子了。

欧拉函数

题目描述

给定n个正整数$a_i$,请你求出每个数的欧拉函数。

输入格式

第一行包含整数n。

接下来n行，每行包含一个正整数$a_i$。

输出格式

输出共n行，每行输出一个正整数$a_i$的欧拉函数。

数据范围

$1 \le n \le 100$,

$1 \le a_i,b_i \le 2*10^9$

输入样例

3
3
6
8

输出样例

2
2
4

实现代码

#include <iostream>

using namespace std;

int main() {
    int n;
    cin >> n;
    while (n--) {
        int x;
        cin >> x;

        int res = x;
        for (int i = 2; i <= x/i; ++i) {
            if (x%i == 0) {
                res = res / i * (i-1);
                while (x%i == 0) x /= i;
            }
        }
        if (x > 1) res = res / x * (x-1);

        cout << res << endl;
    }

    return 0;
}

就是根据定义公式实现代码即可。时间复杂度是$O(\sqrt n)$,主要是在分解质因数上。

筛法求欧拉函数

题目描述

给定一个正整数n,求1～n中每个数的欧拉函数之和。

输入格式

输入一行，包含整数n。

输出格式

共一行，包含一个整数，表示1～n中每个数的欧拉函数之和。

数据范围

$1 \le n \le 10^6$

输入样例

6

输出样例

12

实现代码

#include <iostream>

using namespace std;

typedef long long ll;

const int N = 1e6+10;

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

ll get_eulers(int n) {
    ll res = 0;
    phi[1] = 1;

    for (int i = 2; i <= n; ++i) {
        if (!st[i]) {
            primes[cnt++]  = i;
            phi[i] = i-1; // 如果一个数为质数，那么它与前面的数均互质，因此它的欧拉函数为自身减一
        }
        for (int j = 0; primes[j] <= n/i; ++j) {
            st[primes[j] * i] = true;
            if (i%primes[j] == 0) {
                phi[primes[j]*i] = phi[i] * primes[j];
                break;
            }
            phi[primes[j]*i] = phi[i] * (primes[j]-1);
        }
    }

    for (int i = 1; i <= n; ++i) res += phi[i];

    return res;
}

int main() {
    int n;
    cin >> n;

    cout << get_eulers(n) << endl;

    return 0;
}

其实就是在筛质数的代码基础上进行改动。因为在遍历到一个质数时，这个质数的欧拉函数就是自身减一，即质数的欧拉函数是确定的。

之后再更新i*primes[j]的欧拉函数（这样是求非质数的欧拉函数），此时有两种情况：

1. primes[j]是i的质因数

   由欧拉函数定义，必有：

   $\varPhi(i)=i(1-\frac{1}{p_1})\cdots(1-\frac{1}{p_j})\cdots$。其中 $p_j$ 为primes[j].

   若i可以用基本算数定理表达为:$p_1^{\alpha_1}p_2^{\alpha_2}\cdotsp_j^{\alpha_j}*\cdots$

   则primes[j]*i可表达为 $p_1^{\alpha_1}p_2^{\alpha_2}\cdotsp_j^{\alpha_j+1}*\cdots$​

   由欧拉函数定义可知，欧拉函数只与底数有关，与指数无关，因此 $\varPhi(i)$ 和 $\varPhi(pji)$ 有相同的部分：$(1-\frac{1}{p_1})\cdots(1-\frac{1}{p_j})\cdots$​

   因此:

   $\varPhi(p_ji) = (p_ji)(1-\frac{1}{p_1})\cdots(1-\frac{1}{p_j})\cdots = p_j[i(1-\frac{1}{p_1})\cdots(1-\frac{1}{p_j})\cdots] = p_j\varPhi(i)$

2. primes[j]不是i的质因数

   由欧拉函数定义，必有：

   $\varPhi(i)=i(1-\frac{1}{p_1})\cdots*(1-\frac{1}{p_k})$ 。其中一定不包含$p_j$.

   若i可以用基本算数定理表达为:$p_1^{\alpha_1}p_2^{\alpha_2}\cdotsp_k^{\alpha_k}$，注意，因为pj不是i的质因数，因此一定不包含底数为primes[j]的项

   则primes[j]*i可表达为 $p_1^{\alpha_1}p_2^{\alpha_2}\cdotsp_k^{\alpha_k}*p_j^{\alpha_j}$

   可以发现primes[j]*i相比于i只多了一个质因数（primes[j]）。因此由欧拉函数定义：$\varPhi(p_ji)=(p_ji)(1-\frac{1}{p_1})\cdots(1-\frac{1}{p_k})(1-\frac{1}{p_j})=p_j\varPhi(i)(1-\frac{1}{p_j})=\varPhi(i)*(p_j-1)$

求欧拉函数作用

欧拉定理：若 $a$ 与 $n$ 互质，则 $a^{\varPhi(n)} \mod n = 1$

举例：

a=5，n=6。5和6互质（满足条件）。$\varPhi(6)=2$

$5^2=25，25 \mod 6 = 1$

快速幂

用$O(\log k)$的时间复杂度，快速求出$a^k \mod p$。

如果用传统暴力做法：

int res = 1;
for (int i = 1; i <= k; ++i) res *= a%p;

这样的时间复杂度为 $O(k)$

快速幂算法的核心思想就是：重复利用 $x^y \mod p$，其中 $y$ 为 $2^0、2^1、\cdots、2^{\log k}$

例如：$2^6=64$，如果暴力来做，就是6个2相乘；用快速幂就可以重复利用 $2^{2^1}=4$ ​，即3个4相乘（4由2个2相乘），即最终进行了5次乘法晕眩，少了一次（数量大的时候少的次数更多）

原理其实就是把k拆分由 $2^0、2^1、\cdots、2^{\log k}$ ​中的若干个数相加。即把k用二进制表示，哪些位是1，就加上2的这么多次方。

举例：$4^5 \mod 10$

1. 先预处理出$x^y \mod p，其中y=2^0、2^1、\cdots、2^{\log k}$​

   这里的$\log k$其实就是$\sqrt k$，因此就是$\sqrt 5 \approx 2$​

   $4^{2^0} \mod 10 = 4^1 \mod 10 = 4$

   $4^{2^1} \mod 10 = (4^{2^0} \mod 10)^2 \mod 10=16 \mod 10 =6$

   $4^{2^2} \mod 10 = (4^{2^1} \mod 10)^2 \mod 10 = 36 \mod 10 = 6$

2. 再根据表示进行相乘

   5的二进制表示为101，即由$2^0+2^2$组成，$4^{2^0} \mod 10$和$4^{2^2} \mod 10$的值可查上一步表

   所以$4^5 \mod 10=4^{2^0}4^{2^2} \mod 10=46\mod 10 = 4$

快速幂

题目描述

给定n组$a_i,b_i,p_i$，对于每组数据，求出$a_i^{b_i} \mod p_i$的值。

输入格式

第一行包含整数n。

接下来n行，每行包含三个整数$a_i,b_i,p_i$

输出格式

对于每组数据，输出一个结果，表示$a_i^{b_i} \mod p_i$的值。

每个结果占一行。

数据范围

$1 \le n \le 10^6$​,

$1 \le a_i, b_i, p_i \le2*10^9$

输入样例

2
3 2 5
4 3 9

输出样例

4
1

实现代码

#include <iostream>

using namespace std;

typedef long long ll;

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

int main() {
    int n;
    cin >> n;
    while (n--) {
        int a, k, p;
        cin >> a >> k >> p;
        cout << qmi(a, k, p) << endl;
    }

    return 0;
}

其实这个代码看起来就只是求k的二进制表达。

每一次循环都有a*a%p，即每一次循环结束后，a都变为了$a^2$。第一次循环中的$a=a^{2^0}$；第一次循环结束后，第二次循环中$a=a^{2^1}$；第二次循环结束后，第三次循环中的$a=(a^2)^2=a^{2^2}$，即第i次中的$a=a^{2^{i-1}}$。

每次循环都会判断k的末位，然后判断结束后将k右移一位，相当于第i次循环就是判断k的第i-1位（二进制的数法一般都是从右往左数，且实际的第一位为$2^0$），如果k的第i-1位为1了，可以直接使用a，此时的a表示的就是预处理后的$a^{2^{i-1}}$

快速幂求逆元

题目描述

给定n组$a_i,p_i$，其中$p_i$是质数，求出$a_i \mod p_i$的乘法逆元，若逆元不存在则输出impossible。

乘法逆元定义：

• 若整数$b,m$互质，并且$b|a$，则存在一个整数$x$，使得$a/b \mod m=a*x \mod m$，则称$x$为$b$的模$m$乘法逆元，记为$b^{-1}(mod\ m)$。
• $b$存在乘法逆元的充要条件是$b$与模数$m$互质。当模数$m$为质数时，$b^{m-2}$即为b的乘法逆元。

输入格式

第一行包含整数n。

接下来n行，每行包含两个整数$a_i,p_i$，数据保证$p_i$是质数。

输出格式

共输出n行，每组数据输出一个结果，每个结果占一行。

若$a_i$模$p_i$的乘法逆元存在，则输出一个整数，表示逆元，否则输出impossible。

数据范围

$1 \le n \le 10^5$​,

$1 \le a_i, p_i \le2*10^9$

输入样例

3
4 3
8 5
6 3

输出样例

1
2
impossible

实现代码

其实从题目中的表述可以看出，给出$a_i和p_i$，求$a_i^{p_i-2} \mod p_i$。于是写出如下代码：

#include <iostream>

using namespace std;

typedef long long ll;

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

int main() {
    int n;
    cin >> n;
    while (n--) {
        int a, p;
        cin >> a >> p;
        ll res = qmi(a, p-2, p);
        if (a % p) cout << res << endl;
        else cout << "impossible\n";
    }
}

主要就是改了调用函数的地方。还有就是探讨什么时候无逆元（其实就是在a%p==0时无逆元）

扩展欧几里得算法

「裴蜀定理」

任意一对正整数$a,b$，那么一定存在非零整数$x,y$使得$ax+by=gcd(a, b)$​

证明：

设$gcd(a, b) = d,a=kd,b=md$

因此$ax+by=kxd+myd=(kx+my)*d$，即$ax+by$一定是$a,b$的最小公约数的倍数

接下来就是找一对系数，来满足这样的条件，就需要使用扩展欧几里得算法了。

扩展欧几里得算法

题目描述

给定n对正整数$a_i,b_i$，对于每对数，求出一组$x_i,y_i$，使其满足$a_ix_i+b_iy_i=gcd(a_i， b_i)$。

输入格式

第一行包含整数n。

接下来n行，每行包含两个整数$a_i,b_i$

输出格式

共输出n行，对于每组数据，求出一组满足条件的$x_i,y_i$。每个结果占一行。

本题答案不唯一，输出任意满足条件的$x_i,y_i$均可

数据范围

$1 \le n \le 10^5$​,

$1 \le a_i, b_i \le2*10^9$

输入样例

2
4 6
8 18

输出样例

-1 1
-2 1

实现代码

#include <iostream>

using namespace std;

int exgcd(int a, int b,int& x, int& y) {
    if (!b) {
        x= 1, y = 0;
        return a;
    }

    int d = exgcd(b, a%b, y, x);
    y -= a/b*x;

    return d;
}

int main() {
    int n;
    cin >> n;
    while (n--) {
        int a, b, x, y;
        cin >> a >> b;

        exgcd(a, b, x, y);

        cout << x << " " << y << endl;
    }
    
    return 0;
}

「欧几里得算法」：

int gcd(int a, int b) {
    return b ? gcd(b, a%b) : a;
}

其中证明过了，$gcd(a, b) = gcd(b, a\;mod\ b)$​​

因此拆开写后，递归函数就有两种情况：

1. b==0

   此时$gcd(a,b)=a$

   当且仅当$x=1,y=0$时，$ax+by=a$​

2. b!=0

   int d= exgce(b, a%b, y, x)此时递归函数结束时，一定存在$by+(a\;mod\ b)x=gcd(a, b) = d$​

   $\because a\;mod\ b = a-\lfloor\frac{a}{b}\rfloor*b$​

   $\therefore by+(a-\lfloor\frac{a}{b}\rfloor*b)x=d$​

   $\therefore ax+b(y-\lfloor\frac{a}{b}\rfloor x) = d$​

   所以此时$x$的值不修改，$y-=\lfloor\frac{a}{b}\rfloor x$

线性同余方程

题目描述

给定n对正整数$a_i,b_i,m_i$，对于每组数求出一个$x_i$，使其满足$a_i*x_i\equiv b_i(mod\ m_i)$，如果无解则输出impossible。

输入格式

第一行包含整数n。

接下来n行，每行包含一组数据$a_i,b_i,m_i$

输出格式

共输出n行，每组数据输出一个整数表示一个满足条件的$x_i$，如果无解则输出impossible。

每组数据结果占一行，结果可能不唯一，输出任意一个满足条件的结果均可。

输出答案必须在int范围之内。

数据范围

$1 \le n \le 10^5$​,

$1 \le a_i, b_i, m_i \le2*10^9$

输入样例

2
2 3 6
4 3 5

输出样例

impossible
-3

实现代码

$ax \equiv b(mod\ m) \Leftrightarrow (\exists y \in Z) s.t. ax = my+b$​

即：$ax-my=b$,令$y’=-y$，则$ax+my’=b$，这就是扩展欧几里得。

因此有解的充分必要条件为$gcd(a, m) \mid b$​​

如果直接用exgcd(a, m, x, y)，求出的x满足$ax+my=gcd(a, m)$，但是b不一定等于$gcd(a,m)$，但是有解（即$b=k*gcd(a, m)$）。因此$akx+mky=d$,$x’=kx=\frac{b}{gcd(a,m)}x$.

#include <iostream>

using namespace std;

int exgcd(int a, int b,int& x, int& y) {
    if (!b) {
        x= 1, y = 0;
        return a;
    }

    int d = exgcd(b, a%b, y, x);
    y -= a/b*x;

    return d;
}

int main() {
    int n;
    cin >> n;
    while (n--) {
        int a, b, m, x, y;
        cin >> a >> b >> m;

        int d = exgcd(a, m, x, y);

        if (b % d) cout << "impossible\n";
        else cout << (long long)x*(b/d)%m << endl; 
    }

    return 0;
}

中国剩余定理

中国剩余定理，也叫孙子定理，之所以叫这个名字，是因为《孙子算经》中有这样一个问题：

有物不知其数，三三数之剩二，五五数之剩三，七七数之剩二。问物几何？

这个被叫做“物不知数”的问题本质上是解下面的同余方程组：

后来的数学家在研究中发现，这一方程组有解的一个充分条件是 两两互质，并用构造法给出了这种情况下方程的通解。

即：

$m_1,m_2,\cdots,m_k$两两互质，对于方程组：

$\begin{cases} x \equiv a_1(mod\ m_1) \ x \equiv a_2(mod\ m_2) \ \qquad\quad \vdots \ x \equiv a_k(mod\ m_k) \end{cases}$​

令$M=m_1m_2\cdots*m_k$，$M_i=\frac{M}{m_i}$，$M_i^{-1}$表示$M_i$模$m_i$的逆。

方程组的通解为：$x=a_1M_1M_1^{-1} + a_2M_2M_2^{-1} + \cdots + a_kM_kM_k^{-1}$​

表达整数的奇怪方式

题目描述

给定$2n$个正整数$a_1,a_2,\cdots,a_n$和$m_1,m_2,\cdots,m_n$，求一个最小的整数$x$，满足$\forall i \in [1,n], x \equiv m_i(mod\ a_i)$。

输入格式

第1行包含整数$n$。

第2到$n$行，每$i+1$行包含两个整数$a_i$和$m_i$，数之间用空格隔开。

输出格式

输出整数$x$，如果$x$不存在，则输出-1.

数据范围

$1 \le a_i \le 2^{31}-1$​,

$0 \le m_i \le a_i$

输入样例

2
8 7
11 9

输出样例

31

实现代码

#include <iostream>

using namespace std;

typedef long long ll;

ll exgcd(ll a, ll b, ll &x, ll &y) {
    if (!b) {
        x = 1, y = 0;
        return a;
    }
    
    ll d = exgcd(b, a%b, y, x);
    y -= a/b * x;
    
    return d;
}

int main() {
    int n;
    cin >> n;
    
    bool flag = true;
    ll a1, m1;
    
    cin >> a1 >> m1;
    
    for (int i = 0; i < n - 1; ++i) {
        ll a2, m2;
        cin >> a2 >> m2;
        
        ll k1, k2;
        ll d = exgcd(a1, a2, k1, k2);
        if ((m2-m1) % d) {
            flag = false;
            break;
        }
        
        k1 *= (m2-m1)/d;
        ll t = a2/d;
        k1 = (k1%t + t) % t; // k1可能为负，因此要如此取模
        
        m1 = a1 * k1 * m1;
        a1 = abs(a1/d*a2);
    }
    if (flag) cout << (m1 % a1 + a1) % a1 << endl;
    else cout << "-1\n";
    
    return 0;
}

高斯消元

可以在$O(n^3)$的时间复杂度，求解$n$元线性方程组。

$n$元线性方程组的解有三种情况：

1. 无解

   初等行列变换后，存在0=b的情况

2. 无穷多组解

   初等行列变换后，存在0=0的情况

3. 唯一解

   其他情况

「实现思路」

先枚举每一列，即一列一列地消除：

1. 找到绝对值最大的一行

2. 将该行换到最上面

3. 将该行第一个数变成1

   这一行所有元素同时除以第一个数，实际应该最后改变第一个数的值

4. 将下面所有行的第c列消成0

消除完成后再从下往上，从右往左，将对角线右上角化为0。最终结果为最后一列对应的值。

高斯消元解线性方程组

题目描述

输入一个包含n个方程n个未知数的线性方程组。

$\begin{cases} a{11}x{11} + a{12}x{12} + \cdots + a{1n}x{1n} = b1 \ a{21}x{21} + a{22}x{22} + \cdots + a{2n}x{2n} = b_2 \ \qquad \qquad\qquad\quad \vdots \ a{n1}x{n1} + a{n2}x{n2} + \cdots + a{nn}x_{nn} = b_n \end{cases}$

方程组中的系数为实数。求解这个方程组。

输入格式

第1行包含整数$n$。

接下来n行，每行包含n+1个实数，表示一个方程的n个系数以及等号右侧的常数

输出格式

如果给定线性方程组存在唯一解，则输出共n行，其中第i喊输出第i个未知数的解，结果保留两位小数。

如果给定线性方程组存在无效解，则输出“Infinite group solutions”。

如果给定线性方程组无解，则输出“No solutions”

数据范围

$1 \le n \le 100$，

所有输入系数以及常数均保留两位小数，绝对值均不超过100。

输入样例

3
1.00 2.00 -1.00 -6.00
2.00 1.00 -3.00 -9.00
-1.00 -1.00 2.00 7.00

输出样例

1.00
-2.00
3.00

实现代码

#include <iostream>
#include <cmath>

using namespace std;

const int N = 110;
const double eps = 1e-6;

int n;
double a[N][N];

int gauss() {
    int c, r; // c 表示枚举的列，r表示枚举的行
    for (c = 0, r = 0; c < n; ++c) {
        int t = r;

        // 找到绝对值最大的一行
        for (int i = r; i < n; ++i) if (fabs(a[i][c] > fabs(a[t][c]))) t = i;
        if (fabs(a[t][c]) < eps) continue; // 如果这一行最大的数的绝对值为0，直接枚举下一行

        // 将该行换到当前遍历行的最上一行，即绝对值最大的行t 与 当前行r进行交换
        for (int i = c; i <= n; ++i) swap(a[t][i], a[r][i]);
        
        // 将该行第一个数变成1
        for (int i = n; i >= c; i--) a[r][i] /= a[r][c];
        
        // 将下面所有行的第c列消成0
        for (int i = r + 1; i < n; ++i)
            if (fabs(a[i][c]) > eps) // 当前数不为0才有意义
                for (int j = n; j >= c; --j)
                    a[i][j] -= a[r][j] * a[i][c];

        r++;
    }

    if (r < n) { // 不是唯一解
        for (int i = r; i < n; ++i)
            if (fabs(a[i][n] > eps)) return 2; // 无解
        return 1; // 有无穷多组解
    }

    // 有解，求出解，消除右上角
    for (int i = n-1; i >= 0; --i)
        for (int j = i+1; j < n; ++j)
            a[i][n] -= a[i][j] * a[j][n];

    return 0; // 唯一解
}

int main() {
    cin >> n;
    for (int i = 0; i < n; ++i)
        for (int j = 0;j < n + 1; ++j)
            cin >> a[i][j];

    int t = gauss();

    if (t == 0) { // 唯一解
        for (int i = 0; i < n; ++i) cout << a[i][n] << endl;
    } else if (t == 1) {
        puts("Infinite group solutions");
    } else {
        puts("No solutions");
    }

    return 0;
}

「每一步矩阵的变化情况」：

处理第0行第0列，第一步结束（无变化）
1.00 2.00 -1.00 -6.00
2.00 1.00 -3.00 -9.00
-1.00 -1.00 2.00 7.00

处理第0行第0列，第二步结束（第0行与第1行交换）
2.00 1.00 -3.00 -9.00
1.00 2.00 -1.00 -6.00
-1.00 -1.00 2.00 7.00

处理第0行第0列，第三步结束（第0行整体除a[0][0]）
1.00 0.50 -1.50 -4.50
1.00 2.00 -1.00 -6.00
-1.00 -1.00 2.00 7.00

处理第0行第0列，第四步结束（1～2行的第i行第j列减去a[i][0]*a[0][j]）
1.00 0.50 -1.50 -4.50
0.00 1.50 0.50 -1.50
0.00 -0.50 0.50 2.50

处理第1行第1列，第一步结束（无变化）
1.00 0.50 -1.50 -4.50
0.00 1.50 0.50 -1.50
0.00 -0.50 0.50 2.50

处理第1行第1列，第二步结束（无变化，因为剩下两行里第1行第1列的值本来就是最大的，不需要交换）
1.00 0.50 -1.50 -4.50
0.00 1.50 0.50 -1.50
0.00 -0.50 0.50 2.50

处理第1行第1列，第三步结束（第1行整体除a[1][1]）
1.00 0.50 -1.50 -4.50
0.00 1.00 0.33 -1.00
0.00 -0.50 0.50 2.50

处理第1行第1列，第四步结束（2～2行的第i行第j列减去a[i][1]*a[1][j]））
1.00 0.50 -1.50 -4.50
0.00 1.00 0.33 -1.00
0.00 0.00 0.67 2.00

处理第2行第2列，第一步结束（无变化）
1.00 0.50 -1.50 -4.50
0.00 1.00 0.33 -1.00
0.00 0.00 0.67 2.00

处理第2行第2列，第二步结束（无变化，因为就剩下一行了，该行的第2列一定是最大的，不需要交换）
1.00 0.50 -1.50 -4.50
0.00 1.00 0.33 -1.00
0.00 0.00 0.67 2.00

处理第2行第2列，第三步结束（第2行整体除a[2][2]）
1.00 0.50 -1.50 -4.50
0.00 1.00 0.33 -1.00
0.00 0.00 1.00 3.00

处理第2行第2列，第四步结束（无变化，因为没有其他行的第2列需要被消为0）
1.00 0.50 -1.50 -4.50
0.00 1.00 0.33 -1.00
0.00 0.00 1.00 3.00

最终判断有唯一解时，将右上角消为0后（这里不用实际变为0，只是修改最后一列值即可）
1.00 0.50 -1.50 1.00
0.00 1.00 0.33 -2.00
0.00 0.00 1.00 3.00

高斯消元解异或线性方程组

题目描述

输入一个包含$n$个方程$n$个未知数的异或线性方程组。

方程组中的系数和常数为0或1，每个未知数的取值也为0或1。

求解这个方程组。

异或线性方程组示例如下：

M[1][1]x[1] ^ M[1][2]x[2] ^ ... ^ M[1][n]x[n] = B[1]
M[2][1]x[1] ^ M[2][2]x[2] ^ ... ^ M[2][n]x[n] = B[2]
...
M[n][1]x[1] ^ M[n][2]x[2] ^ ... ^ M[n][n]x[n] = B[n]

输入格式

第一行包含整数n。

接下来n行，每行包含$n+1$个整数0或1，表示一个方程的$n$个系数以及等号右侧的常数。

输出格式

如果给定线性方程组存在唯一解，则输出共$n$行，其中第$i$行输出第$i$个未知数的解。

如果给定线性方程组存在无数解，则输出”Infinite group solutions”。

如果给定线性方程组无解，则输出”No solution”。

数据范围

$1 \le n \le 100$​

输入样例

3
1 1 0 1
0 1 1 0
1 0 0 1

输出样例

1
0
0

实现代码

异或运算等价于不进位的加法。所以可以看成线性方程。

1. 消成上三角矩阵
   1. 枚举一列，找到非零行（非零即一，一即最大）
   2. 把非零行交换到最上面一行
   3. 将该列下面的行消为零（这里和线性方程组的消除方式不同）
2. 判断解
   1. 完美上三角——唯一解
   2. 有矛盾（0=1或1=0）——无解
   3. 无矛盾——无穷解

#include <iostream>
#include <algorithm>

using namespace std;

const int N = 110;

int n;
int a[N][N];

int gauss() {
    int r, c;
    for (r = c = 0; c < n; ++c) {
        // 找到非零行
        int t = r; // t表示非零行
        for (int i = r; i < n; ++i)
            if (a[i][c]) {
                t = i;
                break;
            }
        if (!a[t][c]) continue; // 该列全为0，直接下一列，这样会导致少一次或多次r++，由此来判断是否是完美上三角

        // 把该行(t)交换到最上行(r)
        for (int i = c; i <= n; ++i) swap(a[t][i], a[r][i]);

        // 将下面的行消为0
        for (int i = r+1; i < n; ++i)
            if (a[i][c])
                for (int j = c; j <= n; ++j)
                    a[i][j] ^= a[r][j];

        r++;
    }

    if (r < n) { // 有矛盾，不满足完美三角形
        for (int i = r; i < n; ++i)
            if (a[i][n]) // 有 0=1 的情况
                return 2; // 无解
        return 1; // 无穷多组解
    }

    // 有解，算出解
    for (int i = n-1; i >=0 ; --i)
        for (int j = i + 1; j < n; ++j)
            a[i][n] ^= a[i][j] & a[j][n];

    return 0; // 唯一解
}

int main() {
    cin >> n;
    for (int i = 0; i < n; ++i)
        for (int j = 0; j < n + 1; ++j)
            cin >> a[i][j];

    int res = gauss();

    if (res == 0) 
        for (int i = 0; i < n; ++i) cout << a[i][n] << endl;
    else if (res == 1) puts("Infinite group solutions");
    else puts("No solution");

    return 0;
}

组合数

$C_a^b = \frac{a(a-1)\cdots(a-b+1)}{123\cdots*b} = \frac{a!}{b!(a-b)!}$​

• 求组合数I

  10万询问，$1\le b \le a \le 2000$

  递推，$O(N^2)$

• 求组合数II

  10万询问，$1\le b \le a \le 10^5$

  预处理$i!$和$(i!)^{-1}$，$O(N*logN)$

• 求组合数III

  20询问，$1 \le b \le a \le 10^{18}$​，$1 \le p \le 10^5$​

  卢卡斯定理，$O(p \log N \log p)$​

• 求组合数IV

  不能取模，值很大，需要用高精度。$1 \le b \le a \le 5000$​

求组合数I

题目描述

给定$n$组询问，每组询问给定两个整数$a,b$，请你输出$C_a^b\;mod\ (10^9+7)$的值。

输入格式

第一行包含整数n。

接下来n行，每行包含一组$a$和$b$。

输出格式

共输出n行，每行输出一个询问的解。

数据范围

$1 \le n \le 10^5$​,

$1 \le b \le a \le 2000$

输入样例

3
3 1
5 3
2 2

输出样例

3
10
1

实现代码

根据数据范围可以分析出不能用暴力来做，发现a，b的取值较小，所以可以先预处理出来，这里主要利用递推式（定理）：$Ca^b=C{a-1}^b+C_{a-1}^{b-1}$

有了这个递推式后，直接用DP的做法来实现代码(初始状态为$C_i^0\quad i\in$)：

#include <iostream>

using namespace std;

const int N = 2010, mod = 1e9+7;

int c[N][N];

void init() {
    for (int i = 0; i < N; ++i)
        for (int j = 0; j <= i; ++j)
            if (!j) c[i][j] = 1;
        else c[i][j] = (c[i-1][j] + c[i-1][j-1])%mod;
}

int main() {
    init();
    
    int n;
    cin >> n;
    while (n--) {
        int a, b;
        cin >> a >> b;
        cout << c[a][b] << endl;
    }

    return 0;
}

求组合数II

题目描述

给定$n$组询问，每组询问给定两个整数$a,b$，请你输出$C_a^b\;mod\ (10^9+7)$的值。

输入格式

第一行包含整数n。

接下来n行，每行包含一组$a$和$b$。

输出格式

共输出n行，每行输出一个询问的解。

数据范围

$1 \le n \le 10^5$​,

$1 \le b \le a \le 10^5$

输入样例

3
3 1
5 3
2 2

输出样例

3
10
1

实现代码

这一道题的询问和上一道题一样，也就是说时间复杂度不允许使用暴力方法。同时这道题的a，b的取值范围也比上一题更大，因此也不能再用上一题的递推来做了。

根据定义：$C_a^b=\frac{a!}{(a-b)!b!}$

这道题可以预处理阶乘（$feet[i]=i! \% (10^9+7)$）来做。

但是$\frac{a}{b} \; mod\ p \not = \frac{a\;mod p}{b\;mod\ p}$。因此要用逆元来将乘法转换为乘法（$infact[i]=(i!)^{-1} \% (10^9+7)$）。

因此问题转换为：$C_a^b=fact[a]infact[b-a]infact[b]$

这里可以用快速幂求逆元来求逆元。

#include <iostream>

using namespace std;
typedef long long ll;

const int N = 100010, mod = 1e9+7;

int fact[N], infact[N];

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

int main() {
    // 预处理阶乘和逆元的阶乘
    fact[0] = infact[0] = 1;
    for (int i = 1; i < N; ++i) {
        fact[i] = (ll)fact[i-1]*i%mod;
        infact[i] = (ll)infact[i-1] * qmi(i, mod-2, mod) % mod;
    }

    int n;
    cin >> n;
    while (n--) {
        int a, b;
        cin >> a >> b;
        ll ans = (ll)fact[a]*infact[b]%mod*infact[a-b]%mod;
        cout << ans << endl;
    }

    return 0;
}

求组合数III

题目描述

给定$n$组询问，每组询问给定三个整数$a,b,p$，其中$p$是质数，请你输出$C_a^b\;mod\ p$的值。

输入格式

第一行包含整数n。

接下来n行，每行包含一组$a,b,p$。

输出格式

共输出n行，每行输出一个询问的解。

数据范围

$1 \le n \le 20$​,

$1 \le b \le a \le 10^{18}$​

$1 \le p \le 10^5$

输入样例

3
5 3 7
3 1 5
6 4 13

输出样例

3
3
2

实现代码

「卢卡斯定理」：$Ca^b \equiv C{a\;mod\ p}^{b\; mod\ p}*C_{a/p}^{b/p} (mod\ p)$

#include <iostream>

using namespace std;
typedef long long ll;

int p;

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

int C(int a, int b) {
    int res = 1;
    for (int i = 1, j = a; i <= b; ++i, j--) {
        res = (ll) res * j % p;
        res = (ll)res*qmi(i, p-2)%p;
    }
    return res;
}

int lucas(ll a, ll b) {
    if (a < p && b < p) return C(a, b);
    return (ll) C(a%p, b%p) * lucas(a/p, b/p) %p;
}

int main() {
    int n;
    cin >> n;
    while (n--) {
        ll a, b;
        cin >> a >> b >> p;
        cout << lucas(a, b) << endl;
    }

    return 0;
} 

求组合数IV

题目描述

输入$a,b$，求$C_a^b$的值。

注意结果可能很大，需要使用高精度计算。

输入格式

共一行，包含两个整数$a$和$b$

输出格式

共一行，输出$C_a^b$的值。

数据范围$​

$1 \le b \le a \le 5000$

输入样例

5 3

输出样例

10

实现代码

这道题要用高精度来计算，但是数据范围还是比较大，而且涉及到高精度乘法和高精度除法，因此考虑把高精度除法也转换为乘法。这里就可以使用基本算数定理，将分子分母用质数表示出来，最终给质数进行高精度乘法。

即：$C_a^b=\frac{a!}{(a-b)!b!}$，其中$a! = p_1^{\alpha_1}p_2^{\alpha_2}\cdots p_k^{\alpha_k},(a-b)!b! = p_1^{\beta_1}p_2^{\beta_2}\cdots p_m^{\beta_m}$，

则$C_a^b=p_1^{\alpha_1-\beta_1}p_2^{\alpha_2-\beta_2}\cdots p_m^{\alpha_m-\beta_m}\cdots p_k^{\alpha_a-\alpha_k}$​

其中分子或分母没有的质数的指数当然是0，也可以计算。

「分解$a!$的质因数」：

例如讨论它的某一个质数$p$。则$a! = \lfloor\frac{a}{p}\rfloor + \lfloor\frac{a}{p^2}\rfloor+\cdots+\lfloor\frac{a}{p^k}\rfloor$。$\lfloor\frac{a}{p}\rfloor$表示1～a中p的倍数的个数。

因此本题思路为：

1. 筛质数

   将2～5000中的质数筛选出来

2. 分解质因数

   分解分子和分母的质因数

3. 高精度乘法

#include <iostream>
#include <vector>

using namespace std;

const int N = 5010;

int primes[N], cnt;
bool st[N];
int sum[N]; // 每个质数的次数

void get_primes(int n) {
    for (int i = 2; i <= n; ++i) {
        if (!st[i]) primes[cnt++] = i;
        for (int j = 0; primes[j] <= n/i; ++j) {
            st[primes[j]*i] = true;
            if (i%primes[j] == 0) break;
        }
    }
}

// 求n中包含p的个数
int get(int n, int p) {
    int res = 0;
    while(n) {
        res += n/p;
        n /= p;
    }
    return res;
}

vector<int> mul(vector<int> a, int b) {
    vector<int> c;
    int t = 0;
    for (int i = 0; i < a.size(); ++i) {
        t += a[i] * b;
        c.push_back(t%10);
        t /= 10;
    }
    while (t) {
        c.push_back(t%10);
        t /= 10;
    }
    return c;
}

int main() {
    int a, b;
    cin >> a >> b;

    // 筛质数
    get_primes(a);

    // 分解质因数
    for (int i = 0; i < cnt; ++i) 
        sum[i] = get(a, primes[i]) - get(b, primes[i]) - get(a-b, primes[i]);

    // 高精度乘法
    vector<int> res;
    res.push_back(1);
    for (int i = 0; i < cnt; ++i)
        for (int j = 0; j < sum[i]; ++j)
            res = mul(res, primes[i]);

    // 输出结果
    for (int i = res.size()-1; i>=0; --i) cout << res[i];
    
    return 0;
}

满足条件的01序列

题目描述

给定n个0和n个1，它们将按照某种顺序拍成长度为2n的序列，求它们能排列成的所有序列中，能够满足任意前缀序列中0的个数都不少于1的序列有多少个。

输出的答案对$10^9+7$取模。

输入格式

共一行，包含整数n。

输出格式

共一行，包含一个整数，表示答案。

数据范围$​

$1 \le n \le 10^5$

输入样例

3

输出样例

5

• 提示：

  000111

  001101

  001011

  010011

  010101

实现代码

所有方案数为:$C{2n}^n$，不合法的方案数为$C{2n}^{n-1}$

因此，最终方案数为：$C{2n}^n-C{2n}^{n-1} = \frac{C_{2n}^n}{n+1}$

这样的数被称为「卡特兰数」。

#include <iostream>

using namespace std;
typedef long long ll;

const int mod = 1e9+7;

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

int main() {
    int n;
    cin >> n;

    int a = 2*n, b = n;
    int res = 1;

    // 求a！/(a-b)!
    for (int i = a; i > a-b; --i) res = (ll) res *i%mod;
    // 求b!^{-1}
    for (int i = 1; i <= b; ++i) res = (ll) res*qmi(i, mod-2, mod) %mod;

    res = (ll) res*qmi(n+1, mod-2, mod) % mod;

    cout << res << endl;

    return 0;
}

代码中先求$\frac{a!}{(a-b)!}$,然后用快速幂求1～b的逆元，并累乘，求得$b!^{-1}$

乘得结果为$\frac{a!}{(a-b)!b!} = Ca^b = C{2n}^n$

最后再乘$n+1$的逆元，与前面的结果相乘即可。

容斥原理

百度百科定义

在计数时，必须注意没有重复，没有遗漏。为了使重叠部分不被重复计算，人们研究出一种新的计数方法，这种方法的基本思想是：先不考虑重叠的情况，把包含于某内容中的所有对象的数目先计算出来，然后再把计数时重复计算的数目排斥出去，使得计算的结果既无遗漏又无重复，这种计数的方法称为容斥原理。

以三个集合举例：

要求集合元素个数（VN图包含面积），就是$A \cup B \cup C - A\cap B - B\cap C-A \cap C + A \cap B\cap C$​

能被整除的数

题目描述

给定一个整数$n$和$m$个不同的质数$p_1,p_2,\cdots, p_m$。

请你求出$1$~$n$中能被$p_1,p_2,\cdots,p_m$中的至少一个数整出的整数有多少个

输入格式

第一行包含整数$n$和$m$。

第二行包含$m$个质数。

输出格式

共一行，包含一个整数，表示满足条件的整数的个数。

数据范围

$1 \le m \le 16$,

$1\le n,p_i\le10^9$

输入样例

10 2
2 3

输出样例

7

实现代码

根据样例可以如下转化：

$s_1:{2,4,6,8,10};s_2:{3,6,9}$，即用$s_1$表示2的倍数的集合，$s_2$表示3的倍数的集合。最终要求的结果即为$|s_1 \cup s_2| = |s_1|+|s_2|-|s_1 \cap s_2| = 5 + 3 - 1 = 7$

#include <iostream>

using namespace std;

typedef long long ll;

const int N = 20;

int n, m;
int p[N];

int main() {
    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, cnt = 0;
        for (int j = 0; j < m; ++j) {
            if (i >> j & 1) {
                cnt ++;
                if ((ll)t*p[j] > n) {
                    t = -1;
                    break;
                }
                t *= p[j];
            }
        }
        if (t != -1) {
            if (cnt % 2) res += n/t; // 奇数个集合
            else res -= n/t;
        }
    }

    cout << res << endl;

    return 0;
}

简单博弈论

• NIM游戏

  给定$N$堆物品，第$i$堆物品有$A_i$个。两名玩家轮流行动，每次可以任选一堆，取走任意多个物品，可把一堆取关，但不能不取。取走最后一件物品者获胜。两人都采取最优策略。问娴熟是否必胜。

  我们吧这种游戏称为「NIM博弈」。把游戏过程中面临的状态称为「局面」。整局游戏第一个行动的称为「先手」，第二个行动的称为「后手」。若在某一局面下，无论采取何种行动，都会输掉游戏，则称改局游戏结束。

  所谓采取最优策略是指，若在某一局面下存在某种行动，使得行动后对面面临必败局面，则优先采取该行动。同时，这样的局面被称为「必胜」。我们讨论博弈问题一般都只考虑理想情况。

  NIM博弈不存在平局，至于先手必胜和先手必败两种情况。

  定理：NIM博弈先手必胜，当且仅当$A_1 \oplus A_2 \oplus \cdots \oplus A_n \not= 0$。证明见Nim游戏

• 公平组合游戏ICG

  若一个游戏满足：

  1. 由两名玩家交替行动
  2. 在游戏进程的任意时刻，可以执行的合法行为与轮到哪名玩家无关
  3. 不能行动的玩家判负

  则称该游戏为一个公平组合游戏。

  NIM博弈属于公平组合游戏，但常见的棋类游戏，比如围棋，就是不公平组合游戏。因为围棋交战双方智能落黑子和白子，胜负判定也比较复杂，不能满足第二条和第三条

• 有向图游戏

  给定一个有向无环图，图中有一个唯一的起点，在起点上放有一枚棋子。两名玩家交替地把这枚棋子沿有向边进行移动，每次可以移动一步，无法移动者判负，该游戏称为有向图游戏。

  任何一个公平组合游戏都可以转换为有向图游戏。具体方法是：把每个局面看成图中的一个节点，并且从每个局面沿着合法行动能够到达的下一个局面连有向边。

• Mex运算

  设$S$表示一个非负整数。定义$mex(S)$为求出不属于集合$S$的最小非负整数运算，即：

  $mex(S)=min(x)$，$x\in N,x \not\in S$​

• SG函数

  在有向图游戏中，对于每个节点$x$，设从$x$出发共有$k$条有向边，分别到达节点$y_1,y_2,\cdots,y_k$，定义$SG(x)$为$x$的后继节点$y_1,y_2,\cdots,y_k$的SG函数值构成的集合再执行$mex(S)$运算的结果，即：$SG(x)=mex(SG(y_1),SG(y_2),\cdots,SG(y_k))$

  特别地，整个有向图游戏的SG函数值被定义为有向图游戏起点$s$的SG函数值，即$SG(G)=SG(s)$

Nim游戏

NIM游戏属于公平组合游戏。

如果有$n$堆石子，每堆石子的个数分别为：$a_1, a_2, a_3, a_4,\cdots a_n$，若：$a_1 \oplus a_2 \oplus \cdots \oplus a_n \not= 0$ 则先手必胜，否则先手必输。

证明：

操作到最后的时候有（必败状态）： $0 \oplus 0 \oplus \cdots \oplus 0 = 0$。即每次操作的最优策略就是把剩下的石子的异或值变成0。

• 若当前的异或和不为$0$即$a_1\oplus a_2\oplus \cdots \oplus a_n = x \not=0$,则一定可以通过某种方式使他们的异或和变成$0$，也就是说一定可以通过某种方式从某一堆里面拿走若干个石子，让剩下的这些数的异或值变成$0$

  假设$x$的二进制表示中最高的一位是第$k$位，那么就说明，$a_1$ ~ $a_n$中一定至少存在某一个数$a_i$它的二进制第$k$位为$1$，如果不存在，那很显然不会出现$x$的第$k$位为$1$

  那么就一定会有：$a_i \oplus x \lt a_i$，因为$x$第$k$位之前的位全为$0$与$a_i$进行异或的时候数值不变，但是第$k$位异或后一定变成$0$，所以运算结果会更小

  于是接下来，从$a_i$这一堆中拿走$a_i-(a_i \oplus x)$这么多的石子，那么此时$ai$这一堆的石子个数就变成了$a_i \oplus x$

  那么剩下的所有堆的石子进行异或就是$a_1\oplus a_2\oplus \cdots \oplus ( a_i\oplus x) \oplus\cdots \oplus a_n = x \oplus x = 0$

• 若当前$a_1\oplus a_2 \oplus \cdots \oplus a_n=x=0$，那么无论不管怎么拿，最后得到剩下石子堆进行异或的结果一定不是$0$

  假设从$ai$中拿了一些石子变成$ai’$使得$a_1\oplus a_2 \oplus \cdots \oplus a_i’ \oplus a{i+1} \oplus \cdots \oplus a_n = x = 0$

  将$a1\oplus a_2 \oplus \cdots \oplus a_n = x = 0$与$a_1\oplus a_2 \oplus \cdots \oplus a_i’ \oplus a{i+1} \oplus \cdots \oplus a_n = x = 0$的等式两边分别进行异或，可以得到：$a_i \oplus a_i’ = 0$，但是很显然这两个数是不一样的，所以出现矛盾

  所以不管怎么拿最后剩下的石子异或的结果一定不是$0$

所以当先手遇到的异或的值不为$0$的话，一定可以拿走一些石子使得后手面对的局面异或值为$0$，且后手不管怎么拿留给先手的异或值必然不为$0$，而由于游戏一定可以结束，从而到最后的时候先手拿走最后一些石子，留给后手无石子可拿，先手必赢。

那么反过来，如果先手遇到的异或值为0，则留给后手的异或值一定不是0，后手一定可以拿走一些石子把接下来的异或值为0留给先手，在整个游戏过程中先手遇到的都是异或值为0的状态，到最后先手必输。

题目描述

给定$n$堆石子，两位玩家轮流操作，每次操作可以从任意一堆石子中拿走任意数量的石子（可以拿完，但不能不拿），最后无法进行操作的人视为失败。

问如果两人都采用最优策略，先手是否必胜。

输入格式

第一行包含整数$n$。

第二行包含$n$个数字$a1,a_2,\cdots,a_n$，其中第$i$个数字表示第$i$堆石子的数量。

输出格式

如果先手方必胜，则输出”Yes”。

否则，输出”No”。

数据范围

$1 \le n \le 10^5$,

$1\le a_i\le10^9$

输入样例

2
2 3

输出样例

Yes

实现代码

#include <iostream>

using namespace std;

int main() {
    int n, t, x = 0;
    cin >> n;
    for (int i = 0; i < n; ++i) {
        cin >> t;
        x ^= t;
    }

    if (x != 0) cout << "Yes\n";
    else cout << "No\n";

    return 0;
}

台阶-Nim游戏

题目描述

现在，有一个$n$级台阶的楼梯，每级台阶上都有若干个石子，其中第$i$级台阶上有$a_i$个石子（$i \ge 1$）。

两位玩家轮流操作，每次操作可以从任意一级台阶上拿若干个石子，放到下一级台阶中（不能不拿）。

已经拿到地面上的石子不能再拿，最后无法进行操作的人视为失败。

问如果两人都采用最优策略，先手是否必胜。

输入格式

第一行包含整数$n$。

第二行包含$n$个数字$a1,a_2,\cdots,a_n$，其中第$i$个数字表示第$i$堆石子的数量。

输出格式

如果先手方必胜，则输出”Yes”。

否则，输出”No”。

数据范围

$1 \le n \le 10^5$,

$1\le a_i\le10^9$

输入样例

3
2 1 3

输出样例

Yes

实现代码

根据样例：

最终对手输的情况为：0 0 0，则我的最后一次操作为：从1 0 0变为0 0 0，则我需要保证对手的上一次操作为0 1 0(若对手面对0 3 0时，操作为1 2 0，则我操作为0 2 0，对手再操作为1 1 0，我再操作为0 1 0，所以倒数的情况一定为0 1 0)。即操作后的情况是奇数行的异或值要为0.

这里可以推出结论：$a_1 \oplus a_3\oplus a_5 \oplus \cdots \oplus a_n = x \not= 0$（奇数台阶异或不为零）则先手必胜

#include <iostream>

using namespace std;

int main() {
    int n, t, res = 0;
    cin >> n;
    for (int i = 1; i <= n; ++i) {
        cin >> t;
        if (i & 1) res ^= t;
    }

    if (res) puts("Yes");
    else puts("No");

    return 0;
}

集合-Nim游戏

题目描述

给定$n$堆石子以及一个由$k$个不同正整数构成的数字集合$S$。

现有两位玩家轮流操作，每次操作可以冲任意一堆石子中拿取石子，每次拿取的石子数量必须包含于集合$S$，最后无法进行操作的人视为失败。

问如果两人都采用最优策略，先手是否必胜。

输入格式

第一行包含整数$k$，表示数字集合$S$中数字的个数。

第二行包含$k$个数字，其中第$i$个整数表示数字集合$S$中的第$i$个数$s_i$。

第三行包含整数$n$。

第四行包含$n$个整数，其中第$i$个整数表示第$i$堆石子的数量$h_i$。

输出格式

如果先手方必胜，则输出”Yes”。

否则，输出”No”。

数据范围

$1 \le n,k \le 100$,

$1\le s_i,h_i \le 10000$

输入样例

2
2 5
3
2 4 7

输出样例

Yes

实现代码

根据样例，要拿掉2 4 7 三堆石子，因此可以看成三个有向图游戏，即求三个sg函数值（$sg(2),sg(4),sg(7)$）然后根据Nim游戏的结论（初始异或值不为0则必胜），对三个sg函数值进行异或计算，再判断结果。

前面已经讲过sg函数的计算了，这里举例计算：

先看后继节点：

7->5->3->1

5->0

7->2->0

（7取出2颗变为5，7取出5颗变为2，因此7的后继节点为5和2…）

sg(0) = mex($\varnothing$) = 0 (0没有后继节点)

sg(1) = mex($\varnothing$) = 0 (1没有后继节点)

sg(2) = mex({sg(0)}) = mex({0}) = 1

sg(3) = mex({sg(1)}) = mex({0}) = 1

sg(5) = mex({sg(3), sg(0)}) = mex({1, 0}) = 2

sg(7) = mex({sg(5), sg(2)}) = mex({2, 1}) = 0

到这里求出了sg(7)=2，sg(2)=0

4->2->2->0

sg(4)=mex({sg(2)})=mex({1}) = 0

到这里求出了sg(7)=1，sg(2)=0，sg(4)=0

最终异或值为1。

#include <iostream>
#include <cstring>
#include <algorithm>
#include <unordered_set>

using namespace std;

const int N = 110, M = 10010;

int n, m;
int s[N], f[M]; // s为集合元素，f[i]为sg(i)

int sg(int x) { // 记忆化搜索求f[x]
    if (f[x] != -1) return f[x]; // 保证每个状态只被算一次
    unordered_set<int> S; // 这个集合存储的是sg(y1),sg(y2)...sg(yk)
    for (int i = 0; i < m; ++i) {
        int sum = s[i];
        if (x >= sum) S.insert(sg(x-sum));
    }

    // 从0开始找，第一个不在集合中的数就是mex(S)，即sg(x)
    for (int i = 0; ; ++i) if (!S.count(i)) return f[x] = i;
}

int main() {
    cin >> m;
    for (int i = 0; i < m; ++i) cin >> s[i];

    cin >> n;
    memset(f, -1, sizeof(f));

    int res = 0;
    for (int i = 0; i < n; ++i) {
        int x;
        cin >> x;
        res ^= sg(x);
    }

    if (res) cout << "Yes\n";
    else cout << "No\n";
    
    return 0;
}

拆分-Nim游戏

题目描述

给定$n$堆石子，两位玩家轮流操作，每次操作可以取走其中的一堆石子，然后放入两堆规模更小的石子（新堆规模可以为0，且两个新堆堆石子总数可以大于取走的那堆石子数），最后无法进行操作的人视为失败。

问如果两人都采用最优策略，先手是否必胜。

输入格式

第一行包含整数$n$。

第二行包含$n$个数字$a1,a_2,\cdots,a_n$，其中第$i$个数字表示第$i$堆石子的数量。

输出格式

如果先手方必胜，则输出”Yes”。

否则，输出”No”。

数据范围

$1 \le n,a_i \le 100$

输入样例

2
2 3

输出样例

Yes

实现代码

假设取$a_i$，可以分为两堆$b_i, b_j$，即$sg(a_i) = sg(b_i,b_j) = sg(b_i) \oplus sg(b_j)$

#include <iostream>
#include <cstring>
#include <algorithm>
#include <unordered_set>

using namespace std;

const int N = 110;

int f[N];

int sg(int x) {
    if (f[x] != -1) return f[x];

    unordered_set<int> S;
    for (int i = 0; i < x; ++i) {
        for (int j = 0; j <= i; ++j) {
            S.insert(sg(i) ^ sg(j));
        }
    }

    for (int i = 0; ; ++i)
        if (!S.count(i)) 
            return f[x] = i;
}

int main() {
    int n;
    cin >> n;
    memset(f, -1, sizeof(f));

    int res = 0;
    for (int i = 0; i < n; ++i) {
        int x;
        cin >> x;
        res ^= sg(x);
    }

    if (res) cout << "Yes\n";
    else cout << "No\n";

    return 0;
}