写在开头
- 质数定理：1~n中有n/ln(n)个质数
- 一个数n至多存在一个大于sqrt(n)的质因子
- int范围内，拥有最多约数个数的那个数拥有大约1650个约数，也就是说一个数的约数个数并不是很多
- d|n读作“d能整除n”，意思是n%d==0
- (a,b)读作“a和b的最大公约数”、
- a≡b(mod m)读作“a与b模m同余”，意思是a%m == b%m，a和b除以同一个正整数m后所得余数相同
- a mod b = a - ⌊a/b⌋*b
- d|n ==> (n/d)|n
- d|a 且 d|b ==> d|(a+b)
- d|a 且 d|b ==> d|(a*x+b*y)
- (a,b) = (b,a mod b) // 应用：欧几里得算法
- 互质：两个数的公约数只有1
- 若a与n互质，则a^phi(n)≡1(mod n)
- 若a与c互质，b与c互质，则a*b与c互质
- 因为(a/b) mod p != (a mod p)/(b mod p),所以需要利用逆元把除法运算换成乘法运算
- Iverson bracket 艾弗森括号，是一种用方括号记号，如果方括号内的条件满足则为1，不满足则为0。如[n=0]表示若n=0则取值为1，若n≠0则取值为0
- 写数论题目之前先算时间复杂度，时间复杂度通过了才可以做
质数
判断质数
模板：
bool is_prime(int x)
{
if (x < 2) return false;
for (int i = 2; i <= x / i; i ++ ) // 注意循环条件
if (x % i == 0)
return false;
return true;
}
说明：
- 试除法
- 判断条件的选择
i < n和i <= n / 2时间复杂度都是$O(n)$，过高i * i <= n虽然时间复杂度是$O(n^{\frac{1}{2}})$，但i * i可能会溢出- 因此最好的判别条件是
i <= n / i，时间复杂度是$O(n^{\frac{1}{2}})$
分解质因数
模板：
void divide(int x)
{
for (int i = 2; i <= x / i; i ++ ) // 注意循环条件
if (x % i == 0)
{
int cnt = 0;
while (x % i == 0) x /= i, cnt ++ ;
cout << i << ' ' << cnt << endl;
}
if (x > 1) cout << x << ' ' << 1 << endl; // 至多存在一个大于sqrt(n)的质因子
}
说明：
- 试除法
- 数学知识：$n$最多包含一个大于$\sqrt{n}$的质因子。例如$6=2\times 3$，$\sqrt{6}=2.44949$，存在一个$>\sqrt{6}$的质因子$3$
- 判别质数用
i <= n / i条件 - 最好时间复杂度$O(\text{log}n)$，最坏时间复杂度$O(n^{\frac{1}{2}})$
筛质数
朴素筛法
模板：
int primes[N], cnt; // primes[]存储所有素数（静态链表）
bool st[N]; // st[x]存储x是否被筛掉
void get_primes(){
for(int i = 2; i <= n; i++){
if (!st[i]) {
primes[cnt++]=i; //这里找到的一定是质数，把质数存起来
}
for(int j = i; j <= n; j += i) //不管是合数还是质数，都用来筛掉后面它的倍数
st[j]=true;
}
}
说明：
- 思想：素数的倍数一定不是素数
- 缺点：没必要用合数剔除
- 时间复杂度是$O(n\text{log}n)$
埃氏筛法
模板：
int primes[N], cnt; // primes[]存储所有素数
bool st[N]; // st[x]存储x是否被筛掉
void get_primes(int n)
{
for (int i = 2; i <= n; i ++ )
{
if (!st[i]) {
primes[cnt++] = i; // i是当前可用的最小质数，保存到primes中
for (int j = i + i; j <= n; j += i)
st[j] = true; // 素数的倍数一定不是素数
}
}
}
说明：
- 把朴素筛法的循环放到了if里面，时间复杂度是$O(n\text{loglog}n)$，错误近似可以看成$O(n)$
- 核心思想：只用质数剔除
- 特点：每次发现的第1个非标记数，一定是质数
线性筛法
模板：
int primes[N], cnt; // primes[]存储所有素数
bool st[N]; // st[x]存储x是否被筛掉
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; // 用最小质因子去筛合数primes[j] * i
if (i % primes[j] == 0) break; // 若prime[j]是i的最小质因子，则prime[j+1] * i的最小质因子依旧是prime[j]
}
}
}
说明：
- 核心思想：每个合数必有一个最小素因子，用这个因子筛掉合数
- 理解（见算法基础课四（一）00:34:00之前
- 当
i % primes[j] != 0时，说明此时遍历到的primes[j]不是i的质因子，那么只可能是此时的primes[j] < i的最小质因子，所以primes[j] * i的最小质因子就是primes[j] - 当有
i % primes[j] == 0时，说明i的最小质因子是primes[j]，因此primes[j] * i的最小质因子也就应该是prime[j]，之后接着用st[primes[j+1] * i] = true去筛合数时，就不是用最小质因子去更新了，因为i有最小质因子primes[j] < primes[j+1]，此时的primes[j+1]不是primes[j+1] * i的最小质因子，此时就应该退出循环，避免之后重复进行筛选
- 当
- 时间复杂度是$O(n)$
约数
求所有约数
模板：
vector<int> get_divisors(int x)
{
vector<int> res;
for (int i = 1; i <= x / i; i ++ )
if (x % i == 0)
{
res.push_back(i);
if (i != x / i) res.push_back(x / i);
}
sort(res.begin(), res.end());
return res;
}
说明：
- 试除法求所有约数
- 注意区别因数分解
- 遍历范围是
1~sqrt(x)，每次找到约数时，把自身和另一个同时放入 - 最后需要排序，但实际上消耗的时间不多，
int最大值才有大约1500个约数
约数个数定理
$n = p_1^{c_1}\times p_2^{c_2} *\times … \times p_k^{c_k}$的约数个数等于
$$ (c_1+1)(c_2+1)…(c_k+1) $$
模板：
unordered_map<int, int> primes; // 用哈希表保存质数的指数
// 质数分解
for (int i = 2; i <= x / i; i++)
while(x % i == 0) {
x /= i;
primes[i]++;
}
if (x > 1) primes[x]++;
// 约数个数定理
LL res = 1;
for (auto elem : primes) res = res * (elem.second + 1);
说明：
- 数学原理：
- 一个整数能唯一展开成若干个质数乘积的形式$n = p_1^{c_1}\times p_2^{c_2} \times … \times p_k^{c_k}$
- 整数相乘等价于对应质指数相加$n_1 \times n_2= (p_1^{a_1}\times p_2^{a_2} \times … \times p_k^{a_k}) \times (p_1^{b_1}\times p_2^{b_2} \times … \times p_k^{b_k})= p_1^{a_1+b_1}\times p_2^{a_2+b_2} \times … \times p_k^{a_k+b_k}$
- 例子：$12=2^2\times3^1$，$12$的约数有$1, 2, 3, 4, 6, 12$共$6$个，根据公式计算同样是$(2+1)\times(1+1)=6$个
约数之和定理
$$ sum=\prod_{i=1}^k{\sum_{j=0}^{c_i}{p_{i}^{j}}}=(p_1^0 + p_1^1 + … + p_1^{c_1})\times…\times(p_k^0 + p_k^1 + … + p_k^{c_k}) $$
模板：
LL res = 1;
for (auto elem : primes) {
int p = elem.first, a = elem.second;
LL sum = 1;
while(a--) sum = sum * p + 1;
res *= sum;
}
说明：
-
例子：$12=2^2\times3^1$，$12$的约数有$1, 2, 3, 4, 6, 12$，约数之和为$28$，根据公式计算同样是$(2^0+2^1+2^2)\times(3^0+3^1)=28$个
-
$p^0+p^1+p^2+…+p^n=…((1 \times p + 1)\times p + 1)…$，其中里边迭代
n次。因此可直接用结果进行迭代计算，而不用一个变量存储中间值
最大公约数
模板：
int gcd(int a, int b)
{
return b ? gcd(b, a % b) : a;
}
说明：
- 欧几里得算法
- $\text{gcd}(a, b) = \text{gcd} (b, a\mod b)$
- $\text{gcd} (a, 0) = a$
欧拉函数
$\varphi (n)$表示$1$~$n$中与$n$互质的数的个数。
$$ \varphi \left( n \right) =n\left( 1-\frac{1}{p_1} \right) \left( 1-\frac{1}{p_2} \right) \cdots \left( 1-\frac{1}{p_k} \right) $$
其中$n=p_{1}^{c_1}p_{2}^{c_2}\cdots p_{k}^{c_k}$
求单个欧拉函数
模板：
int phi(int 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);
return res;
}
说明：
- 为避免浮点数除法，编程时采用另一个形式编写代码$\varphi \left( n \right) =n\left( \frac{p_1-1}{p_1} \right) \left( \frac{p_2-1}{p_2} \right) \cdots \left( \frac{p_k-1}{p_k} \right)$
- 可用先除后乘的方式
res = res / n * (n-1)避免结果溢出 - 例子：$\varphi(6)=2$，因为$1$~$6$中只有$1$和$5$与$6$互质
- 可用容斥原理证明：$\varphi \left( n \right) =n-\sum_i{\frac{n}{p_i}+\sum_{i,j}{\frac{n}{p_i p_j}}-\sum_{i,j,k}{\frac{n}{p_i p_j p_k}}}+\cdots$
- 时间复杂度$O(\sqrt{n})$
筛法求多个欧拉函数
模板：
int primes[N], cnt; // primes[]存储所有素数
int euler[N]; // 存储每个数的欧拉函数
bool st[N]; // st[x]存储x是否被筛掉（合数标记）
void get_eulers(int n)
{
euler[1] = 1;
for (int i = 2; i <= n; i ++ )
{
if (!st[i])
{
// 质数
primes[cnt ++ ] = i;
euler[i] = i - 1;
}
for (int j = 0; primes[j] <= n / i; j ++ )
{
int t = primes[j] * i;
st[t] = true;
if (i % primes[j] == 0)
{
euler[t] = euler[i] * primes[j];
break;
}
euler[t] = euler[i] * (primes[j] - 1);
}
}
}
说明：
- 核心思想：$\varphi(n)$与$n$的质指数大小无关
- 当$n$是质数时，$\varphi(n)=n-1$
- 当$n \mod p = 0$时，$\varphi(pn)=p\varphi(n)$，因为质数$p$只影响$pn$的质指数$c_p$
- 当$n \mod p \ne 0$时，$\varphi(pn)=p\varphi(n)(1-\frac{1}{p})=(p-1)\varphi(n)$，因为质数$p$影响了$pn$展开项数，多了一项$p^1$
- 借助线性筛法求欧拉函数，时间复杂度$O(n)$
快速幂
快速计算$a^k \mod p$
模板：
typedef long long LL;
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;
}
说明：
- 实际上是把$k$表示成二进制$b_{\text{log}k}\cdots b_2b_1b_0$，因此$k=b_02^{2^0} + b_12^{2^1}+b_22^{2^2}+ \cdots b_{\text{log}k} 2^{2^{\text{log}k}}$
- 计算$a^k \mod p=a^{b_02^{2^0} + b_12^{2^1}+b_22^{2^2}+ \cdots b_{\text{log}k} 2^{2^{\text{log}k}}} \mod p=(a^{b_0 2^{2^0}}\mod p) \times (a^{b_1 2^{2^1}}\mod p) \times \cdots \times (a^{b_{\text{log}k} 2^{2^{\text{log}k}}}\mod p)$
- 此外迭代满足如下关系$2^{2^{i+1}} \mod p=(2^{2^{i}} \mod p)^2$
- 中间结果可能会溢出，要强制进行类型转换
- 时间复杂度$O(\text{log}k)$
扩展欧几里得算法
模板
// 求x, y，使得ax + by = gcd(a, b)
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 = x'，x = y'
y -= (a/b) * x; // y = x'-(a / b) * y' = y - (a / b) * x
return d;
}
// 详细版
// 求x, y，使得ax + by = gcd(a, b)
int exgcd(int a, int b, int &x, int &y)
{
if (!b)
{
x = 1; y = 0;
return a;
}
int x_new, y_new;
int res = exgcd(b, a % b, x_new, y_new);
x = y_new; // x = y'
y = x_new - (a/b) * x; // y = x' - (a / b) * y'
return res;
}
说明：
- 在计算最大公约数时，顺便求出一组整数解$(x,y)$，使其满足$ax+by=\text{gcd}(a, b)$
- 拓展：方程$ax+by=d$有解的充要条件是$\text{gcd}(a, b) | d$
- 可用于求同余方程$ax\equiv d\left( \mathrm{mod}b \right) \Leftrightarrow \exists y\in \mathrm{Z},s.t. ax+by=d$，但注意，用
ex_gcd(a, b, x, y)求得的x不是最终x，因为此时ex_gcd求出的x是满足的$ax+by=\text{gcd}(a,b)$的，而不是满足$ax+by=d$的，二者相差$d / \text{gcd}(a, b)$倍 - 模数m不是质数时，求逆元只能使用扩展欧几里得算法，如果模数是质数，那么可以利用费马小定理将求逆元转化为快速幂
证明
$$ \begin{aligned} ax+by&=\mathrm{gcd}\left( a,b \right) \\\ &=gcd\left( b,a\mathrm{mod}b \right) \\\ &=bx\prime+a\mathrm{mod}b\,\,y\prime \\\ &=bx\prime+\left( a-\lfloor \frac{a}{b} \rfloor b \right) y\prime \\\ &=bx\prime+ay\prime-\lfloor \frac{a}{b} \rfloor by\prime \\\ &=ay\prime+b\left( x\prime-\lfloor \frac{a}{b} \rfloor b \right) \end{aligned} $$
对比$a$和$b$的系数得：
$$
\begin{cases} x=y\prime\\\ y=x\prime-\lfloor \frac{a}{b} \rfloor y\prime\\\ \end{cases}
$$
因此可根据递归得到的$x\prime$和$y\prime$计算$x$和$y$
中国剩余定理
已知$m_1,m_2\cdots m_k$互质，方程
$$ \left\{ \begin{array}{c} x\equiv a_1\left( \mod m_1 \right)\\\ \begin{array}{l} x\equiv a_2\left( \mod m_2 \right)\\\ \cdots\\\ x\equiv a_k\left( \mod m_k \right)\\\ \end{array}\\\ \end{array} \right. $$
的最小正整数解$x=a_1MM_1^{-1}+a_2MM_2^{-2}+\cdots a_kMM_k^{-k}$，其中$M=m_1m_2\cdots m_k$，$M_i=\frac{M}{m_i}$，$M_i^{-1}$表示$M_i$模$m_i$的逆元，可通过$M_ix\equiv1(\mod m_i)$求出
模板
typedef long long LL;
// 扩展欧几里得算法
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 = x'，x = y'
y -= (a/b) * x; // y = x'-(a / b) * y' = y - (a / b) * x
return d;
}
// 中国剩余定理
LL x = 0, m1, a1;
cin >> m1 >> a1;
for (int i = 0; i < n - 1; i ++ )
{
LL m2, a2;
cin >> m2 >> a2;
LL k1, k2;
LL d = exgcd(m1, -m2, k1, k2); // 求的不是最终解k1和k2
if ((a2 - a1) % d)
{
x = -1;
break; // 无解
}
k1 *= (a2 - a1) / d; // 变换成真正的解
int t = m2 / d; // k1 = k1 + k * m2 / d
k1 = (k1 % t + t ) % t; // 变换成最小的正整数（防止C++对负数模取余）
x = k1 * m1 + a1;
// 把两个方程合并成一个方程
LL m = abs(m1 / d * m2); // 变成正数
a1 = k1 * m1 + a1;
m1 = m;
}
if (x != -1) x = (x % m1 + m1) % m1; // 变成最小正整数（防止C++对负数模取余）
cout << x << endl;
说明
- 本质是求解同余方程组
- 如果直接按公式求解，结果可能会溢出，因此采用迭代的方式求解
- $x\equiv a\left( \mathrm{mod} m \right) \Leftrightarrow x=mk+a,k\in \mathrm{Z}$
- 考虑第1个和第2个方程，可得等式$m_1k_1+a_1=m_2k_2+a_2$，变形得$m_1k_1-m_2k_2=a_2-a_1$
- 令$d=\text{gcd}(m_1, -m_2)$方程有解的充分必要条件是$d | (a_2-a_1)$
- 不定方程的一个通解是$\begin{cases} k_1=k_1+k\frac{m_2}{d}\\ k_2=k_2+k\frac{m_1}{d}\\ \end{cases}$，其中$k$是参数
- 此时$x=k_1m_1+a_1=(k_1+k\frac{m_2}{d})m_1+a_1=k_1m_1+k\frac{m_1m_2}{d}+a_1=k_1m_1+a_1+k\frac{m_1m_2}{d}=x_1+kd_{12}$，这里令$x_1=k_1m_1+a_1$，$d_{12}=\frac{m_1m_2}{d}$
- 若把$x_1$看成$a_1$，$d_{12}$看成$m_1$，则上式变成$x=a_1+km_1$，即$x\equiv a_1 (\mod m_1)$，相当于把两个同余方程合并成了一个
高斯消元
模板
// a[N][N+1]是增广矩阵
int gauss()
{
int c, r;
// 按列遍历，化成行阶梯矩阵
for (c = 0, r = 0; c < n; c ++)
{
// 找到当前列绝对值最大元素所在的行（搜索第r行~最后一行）
int t = r;
for (int i = r + 1; i < n; i ++ )
if (fabs(a[i][c]) > fabs(a[t][c]))
t = i;
if (fabs(a[t][c]) < eps) continue; // 当前列全是0
for (int j = c; j <= n; j ++ ) swap(a[t][j], a[r][j]); // 将绝对值最大的行换到最顶端
for (int j = n; j >= c; j -- ) a[r][j] /= a[r][c]; // 将当前上的首位变成1，注意从后往前遍历
for (int i = r + 1; i < n; i ++ ) // 用当前行将下面所有的列消成0
if (fabs(a[i][c]) > eps)
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; // 出现0!=0，无解
return 1; // 都是0=0，有无穷多组解
}
// 化成单位阵，增广矩阵的扩展部分为方程的解
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; // 有唯一解
}
说明
- 可在$O(n^3)$求解线性方程组
- 找到绝对值最大的一行
- 将该行换到最上面
- 将该行第一个数化成$1$
- 把该列下边的数化成$0$
- 注意部分过程需要从后往前遍历，否则出现写后读问题
组合数
适用情况：
- 递归法：询问次数多但组合数取值范围小
- 乘法逆元：询问次数多，组合数取值范围大，模数为大质数
- 卢卡斯定理：询问次数少，组合数取值范围特别大，模数为质数
- 分解质因数：模数不是质数甚至不取模
递归法（O(n^2)
$$ C_{a}^{b}=C_{a-1}^{b}+C_{a-1}^{b-1} $$
模板
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];
说明
- 意义： $C_{a}^{b}$表示从$a$个苹果中选$b$个的方案数
i=0时，不会执行else部分，因此不会出现数组越界- 代码结构有点像杨辉三角
- 适用于询问次数多，但组合数取值范围小的情况（$10^3$数量级）
预处理逆元（O(nlogn)
假设$m$是一个很大的质数
$(n-1)!$的模$m$逆元是$((n-1)!)^{m-2}$
$n!$的模$m$逆元是$(n!)^{m-2}$
则$\text{infact}(n)=\text{infact}(n-1) \times \frac{(n!)^{m-2}}{((n-1)!)^{m-2}}=\text{infact}(n-1) \times n^{m-2}$
其中$n^{m-2}$可通过快速幂求得
模板
// 快速幂模板
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;
}
// 预处理阶乘的余数和阶乘逆元的余数
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;
}
说明
- 在涉及乘法运算时，注意先强制转换成
long long，然后再取模，否则结果会溢出 - 适用于询问次数多，但组合数取值范围较大的情况（$10^6$数量级）
Lucas定理（O(plognlogp)
若$p$是质数，则对于任意整数 $1 \le m \le n$，有：
$$ C_n^m=C_{n \mod p}^{m \mod p} \times C_{n \div p}^{m \div p} \mod p $$
模板
// 快速幂模板
int qmi(int a, int k)
{
int res = 1;
while (k)
{
if (k & 1) res = (LL)res * a % p;
a = (LL)a * a % p;
k >>= 1;
}
return res;
}
// 通过定理求组合数C(a, b)
int C(int a, int b)
{
int res = 1;
for (int i = 1, j = a; i <= b; i ++, j -- )
{
res = (LL)res * j % p; // 构造分子a * (a - 1) * ... * (a - b + 1)
res = (LL)res * qmi(i, p - 2) % p; // 构造(b!)的逆元(b!)^{p-2}
}
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;
}
说明
- 适用于询问次数少，但组合数取值范围较大的情况（$a, b \le10^{18}$，$p\le 10^5$）
- 两方面优化组合数的求解
- 快速幂加速计算小规模组合数
- 卢卡斯定理降低组合数的规模
- 注意变量类型以及输入输出的类型是
int还是long long
分解质因数法
模板
int primes[N], cnt; // 存储所有质数
int sum[N]; // 存储每个质数的次数
bool st[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;
}
get_primes(a); // 预处理范围内的所有质数
for (int i = 0; i < cnt; i ++ ) // 求每个质因数的次数
{
int p = primes[i];
sum[i] = get(a, p) - get(b, p) - get(a - b, p);
}
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]);
说明
- 本质是质因数分解+高精度乘法
- 适用求真实值，而非取余后的结果
- 主要思想
- 线性筛法找出$2$~$a$内的所有质数
- 依次遍历质数，求出组合数在各个质数的次数（通过
get方法计算） - 高精度乘法计算各个质因子的乘积
get(int n, int p)求n!质因数分解在$p$的次数的数学公式$\lfloor \frac{n}{p} \rfloor +\lfloor \frac{n}{p^2} \rfloor +\lfloor \frac{n}{p^3} \rfloor +\cdots$
卡特兰数
$\text{catalan}(n)=\frac{C_{2n}^n}{n+1}$
容斥原理
$$ \left| \bigcup_{k=1}^n{S_k} \right|=\sum_i{\left| S_i \right|}-\sum_{i,j}{\left| S_i\cap S_j \right|}+\sum_{i,j,k}{\left| S_i\cap S_j\cap S_k \right|}-\cdots $$
$$ C_{n}^{1}+C_{n}^{2}+C_{n}^{3}+\cdots +C_{n}^{n}=2^n-C_{n}^{0}=2^n-1 $$
模板
int res = 0;
// 一共有2^m-1种集合
for (int i = 1; i < 1 << m; i ++ )
{
int t = 1, s = 0; // t表示当前集合的乘积，s表示符号（正负交替）
// 遍历每种集合
for (int j = 0; j < m; j ++ )
if (i >> j & 1)
{
// 元素存在
if ((LL)t * p[j] > n)
{
t = -1; // 构造的质数乘积比n大，舍弃
break; //构造的质数乘积太大意味着不可能出现同时被这些两两互质的质数整除的情况
}
t *= p[j];
s ++ ;
}
if (t != -1)
{
if (s % 2) res += n / t; // 根据s的奇偶性实现正负交替
else res -= n / t;
}
}
说明
- $2^m$可以用
1 << m实现 - m个元素可以构造$2^m-1$种不同的非空集合，集合包含的元素可用$m$位二进制表示，例如
5=101b，表示集合有第1个元素和第3个元素 - 对于每个元素都有$C_{n}^{1}-C_{n}^{2}+C_{n}^{3}-\cdots +\left( -1 \right) ^{k-1}C_{n}^{n}=1$，即每个元素取到的次数都是1，不重复
- $1$~$n$中能被$p$整除的数的个数为$\lfloor \frac{n}{p} \rfloor$
- 遍历的次序与公式不太一样，不是先加$m$个数，再减去$C_m^2$个数。实际上，当集合元素个数是奇数时，符号为正，反之为负。因此可根据集合个数的奇偶性判断符号的正负
-
模板是基于整除个数问题设计的
-
时间复杂度为$O(2^n)$
NIM游戏
公平组合游戏ICG
若一个游戏满足以下条件，则称该游戏为一个公平组合游戏。
- 由两名玩家交替行动；
- 在游戏进程的任意时刻，可以执行的合法行动与轮到哪名玩家无关；
- 不能行动的玩家判负；
NIM博弈属于公平组合游戏，但城建的棋类游戏，比如围棋，就不是公平组合游戏。因为围棋交战双方分别只能落黑子和白子，胜负判定也比较复杂，不满足条件2和条件3。
取石子
给定N堆物品，第$i$堆物品有$A_i$个。两名玩家轮流行动，每次可以任选一堆，取走任意多个物品，可把一堆取光，但不能不取。取走最后一件物品者获胜。两人都采取最优策略，问先手是否必胜。
定义
我们把这种游戏称为NIM博弈。把游戏过程中面临的状态称为局面。整局游戏第一个行动的称为先手，第二个行动的称为后手。若在某一局面下无论采取何种行动，都会输掉游戏，则称该局面必败。
所谓采取最优策略是指，若在某一局面下存在某种行动，使得行动后对面面临必败局面，则优先采取该行动。同时，这样的局面被称为必胜。我们讨论的博弈问题一般都只考虑理想情况，即两人均无失误，都采取最优策略行动时游戏的结果。
NIM博弈不存在平局，只有先手必胜和先手必败两种情况。
必胜状态，先手进行某一个操作，留给后手是一个必败状态时，对于先手来说是一个必胜状态。即先手可以走到某一个必败状态。
必败状态，先手无论如何操作，留给后手都是一个必胜状态时，对于先手来说是一个必败状态。即先手走不到任何一个必败状态。
定理
NIM博弈先手必胜，当且仅当$A_1\oplus A_2\oplus \cdots \oplus A_n\ne 0$
性质
- 操作到最后时，每堆石子数都是0，$0\oplus 0\oplus \cdots \oplus 0=0$
- 在操作过程中，如果 $A_1\oplus A_2\oplus \cdots \oplus A_n=x\ne 0$。那么玩家必然可以通过拿走某一堆若干个石子将异或结果变为0。
- 在操作过程中，如果 $A_1\oplus A_2\oplus \cdots \oplus A_n= 0$，那么无论玩家怎么拿，必然会导致最终异或结果不为$0$
性质证明
- 性质2的证明：不妨设$x$的二进制表示中最高一位$1$在第$k$位，那么在$A_1,A_2,\cdots,A_n$中，必然有一个数$A_i$，它的第$k$为时$1$，且$A_i\oplus x<A_i$，那么从第$i$堆石子中拿走$A_i-A_i\oplus x$个石子，第ii堆石子还剩$A_i-\left( A_i-A_i\oplus x \right) $ ，此时$A_1\oplus A_2\oplus \cdots \oplus A_i\oplus x\oplus \cdots \oplus A_n=x\oplus x=0$
- 性质3的证明：反证法：假设玩家从第$i$堆石子拿走若干个，结果仍是$0$。不妨设还剩下$A_i′$个，因为不能不拿，所以$0≤A_i′<A_i$，且$A_1\oplus A_2\oplus \cdots \oplus A_i′\oplus \cdots \oplus A_n=0$。那么$\left( A_1\oplus A_2\oplus \cdots \oplus A_i\oplus \cdots \oplus A_n \right) \oplus \left( A_1\oplus A_2\oplus \cdots \oplus A_i′\oplus \cdots \oplus A_n \right) =A_i\oplus A_i′=0$，则 $A_i= A_i′$，与假设$0≤A_i′<A_i$矛盾。
分析
- 当先手面对的局面是$A_1\oplus A_2\oplus \cdots \oplus A_n\ne 0$，则先手总有办法让后手局面变成$A_1\oplus A_2\oplus \cdots \oplus A_n= 0$，后手必败，先手必赢
- 当先手面对的局面是$A_1\oplus A_2\oplus \cdots \oplus A_n= 0$，则后手总有办法让先手局面变成$A_1\oplus A_2\oplus \cdots \oplus A_n\ne 0$，先手必败，后手必赢
取石子（升级版）
约束条件
每次可拿石子的个数不是任意的，而是集合$S$中的元素
概念
$\text{mex}$
$\text{mex(S)}$表示不属于集合$S$的最小自然数，即$\text{mex}(S)=\min \left\{ x \mid x\notin S\land x\in \mathrm{N} \right\} $ 。
- 例如当$S={1, 2}$时，$\text{mex}(S)=0$
- 例如当$S={0, 2}$时，$\text{mex}(S)=1$
有向图游戏
给定一个有向无环图，图中有一个唯一的起点，在起点上放有一枚棋子。两名玩家交替地把这枚棋子沿有向边进行移动，每次可以移动一步，无法移动者判负。该游戏被称为有向图游戏。
任何一个公平组合游戏都可以转化为有向图游戏。具体方法是，把每个局面看成图中的一个节点，并且从每个局面向沿着合法行动能够到达的下一个局面连有向边。
$\text{SG}$
在有向图游戏中，对于每个节点$x$，设从$x$出发共有$k$条有向边，分别到达节点$y_1, y_2,\cdots, y_k$，定义$\text{SG}(x)$为$x$的后继节点$y_1, y_2,\cdots, y_k$的$\text{SG}$函数值构成的集合再执行$\text{mex}(S)$运算的结果，即：
$$ \text{SG}(x) = \text{mex}({\text{SG}(y1), \text{SG}(y2),\cdots, \text{SG}(yk)}) $$
特别地，整个有向图游戏$G$的$\text{SG}$函数值被定义为有向图游戏起点$s$的$\text{SG}$函数值，即$\text{SG}(G) = \text{SG}(s)$。
有向图游戏的和
设$G_1, G_2, \cdots, G_m$是$m$个有向图游戏。定义有向图游戏$G$，它的行动规则是任选某个有向图游戏$G_i$，并在$G_i$上行动一步。$G$被称为有向图游戏$G_1, G_2,\cdots, G_m$的和。
有向图游戏的和的SG函数值等于它包含的各个子游戏SG函数值的异或和，即：$$ \mathrm{SG}\left( G \right) =\mathrm{SG}\left( G_1 \right) \oplus \mathrm{SG}\left( G_2 \right) \oplus \cdots \oplus \mathrm{SG}\left( G_m \right) $$
性质
- $\text{SG}(k)有$$k$个后继结点，分别是$0$~$k-1$
- 非$0$可以走向$0$
- $0$只能走向非$0$
性质证明
类似NIM游戏的证明
定理
有向图游戏的某个局面必胜，当且仅当该局面对应节点的SG函数值大于$0$。
有向图游戏的某个局面必败，当且仅当该局面对应节点的SG函数值等于$0$。
模板
int s[N], f[M];
memset(f, -1, sizeof f); // 初始化
// 记忆化搜索（备忘录法）
int sg(int x)
{
if (f[x] != -1) return f[x]; // 已经计算过
// 构造子树
unordered_set<int> S;
for (int i = 0; i < m; i ++ )
{
int sum = s[i];
if (x >= sum) S.insert(sg(x - sum)); // 保存后继结点到S中（递归）
}
// mem(x)
for (int i = 0; ; i ++ )
if (!S.count(i))
return f[x] = i;
}
// 把n堆石子看成n个独立的有向图G，把各个有向图结果做异或即可得到答案
int x, res = 0;
for (int i = 0; i < n; i ++ )
{
cin >> x;
res ^= sg(x);
}
参考资料：y总直播，站内卢盼盼笔记
