数学笔记
数论
线性筛
枚举倍数
int p[N]， f[N]， idx;
void sieve(){
n = read()， f[0] = f[1] = 1;
for(int i = 2;i <= n;i++){
if(!f[i]) p[++idx] = i;
for(int j = 1;j <= idx && i*p[j] <= n;j++){
f[i*p[j]] = 1;
if(i%p[j] == 0) break;
}
}
}
快速幂
\(O(logn)\) 求出 \(a^b\) \(\%m\)
ll qpow(ll x, ll y, ll m){
ll ans = 1, base = x;
while(y){
if(y&1) ans = ans*base%m;
base = base*base%m;
y >>= 1;
}
return ans;
}
欧拉函数
\(\Phi(n)\) 即为从 \(1-n\) 的正整数与 \(n\) 互质的数的个数 通常用 \(phi[n]\) 表示
一些性质：
若 \(n\) 为质数，则 \(phi[n]\) = \(n-1\)
若 \(n = p^k\)（\(p\) 为质数，则 \(phi[n] = p^k - p^{k-1}\)
若 \(n\) 与 \(m\) 互质，则 \(phi[n\times m] = phi[n] \times phi[m]\) （积性函数）
若 \(n\) 为奇数，则 \(phi[2n] = phi[n]\)
设 \(p_i\) 为 \(n\) 的质因数，\(phi[n] = n\times (1-\frac{1}{p_1})\times (1-\frac{1}{p_2})\times (1-\frac{1}{p_{...}})\times (1-\frac{1}{p_n})\)
关于欧拉函数积性函数的性质：
前置：\(p|n\) 指 \(n\) 能整除 \(p\)，即 \(n\%p\) == 0
若 \(p|n\) 且 \(p^2|n\)，则 \(phi[n] = phi[n/p]\times p\)
若 \(p|n\) 且 \(p^2!|n\)，则 \(phi[n] = phi[n/p]\times (p-1)\)
\(\sum_{d|n}phi[d] = n\)
通过分解因数 \(O(\sqrt n)\) 求 \(n\) 的欧拉函数：
int phi(int x){
int res = x;
for(int i = 2;i <= sqrt(x);i++){
if(x%i == 0){
res = res*(i-1)/i;
while(x%i == 0) x /= i;
}
}
if(x > 1) res = res*(x-1)/x;
return res;
}
线性筛 \(O(n)\) 求 \(1-n\) 的欧拉函数：
for(int i = 2;i <= n;i++){
if(!f[i]) p[++idx] = i, phi[i] = i-1;
for(int j = 1;j <= idx && p[j]*i <= n;j++){
f[i*p[j]] = 1;
if(i%p[j] == 0){
phi[i*p[j]] = p[j]*phi[i];
break;
}
else phi[i*p[j]] = (p[j]-1)*phi[i];
}
}
欧几里得算法（gcd）
快速求 gcd，时间复杂度为 \(O(logn)\)
int gcd(int x, int y){
if(y) return gcd(y, x%y);
return x;
}
扩展欧几里得算法（exgcd）
裴蜀定理与 exgcd
裴蜀定理：对于任意整数 a, b 存在 x, y 满足 ax+by = gcd(a, b) x, y 不止一组
推论：若 ax+by = d 有解，则 gcd(a, b)|d。ax+by = 1 有解是 a, b 互质的充要条件
证明：
\(ax_1+by_1 = gcd(a,b)\)
\(bx_2+(a\%b)y_2 = gcd(b, a\%b)\)
又 \(gcd(a,b) = gcd(b,a\%b)\)
移项可得
\(x_1 = y_2\)
\(y_1 = x_2-(a/b)\times y\)
向下递归，当 \(b=0\)，\(gcd(a,b) = a\)，此时 \(x = 1, y = 0\)，返回 \(gcd(a,b)\)，通过递归，求出一组特解（\(x_0,y_0\)）
时间复杂度同普通 gcd，为 \(O(logn)\)
int exgcd(int a, int b, int &x, int &y){
if(b == 0){
x = 1, y = 0;
return a;
}
int r = exgcd(b, a%b, x, y)， t = x;
x = y;
y = t-(a/b)*y;
return r;
}
exgcd 求线性不定方程
对于 \(ax+by = c\)，根据裴蜀定理，若方程存在整数解，需要满足条件：\(gcd(a,b)|c\)。否则，无整数解。
exgcd 求出的一组特解 \(x_0,y_0\)，乘上系数 \(c/gcd(a,b)\)，就从 \(ax+by=gcd(a,b)\) 的解转化到了 \(ax+by=c\) 的解，即 \(x_0\times c/gcd(a,b)\) 和 \(y_0\times c/gcd(a,b)\)。
再加减系数（其实化简后就是加减了 \(a\times b/gcd(a,b)\times t\)），通解公式：
\(x = c/gcd(a,b)\times x_0+b/gcd(a,b)\times t\)
\(y = c/gcd(a,b)\times y_0-a/gcd(a,b)\times t\)
(\(t\in Z\))
int cal(int a, int b, int &x, int &y){
int g = exgcd(a, b, x, y);
if(c%g != 0) return 0;
int k = c/g;
x *= k, y *= k;//此为一组解 所有解为 x*k+b/k*t, y*k+a/k*t (t为整数)
return 1;
}
同余
同余概念
定义：两个整数 \(x,y\)，若他们除以正整数 \(m\) 的余数相等，则称 \(x,y\) 模 \(m\) 同余，记为 \(x\equiv y \pmod m\)
同余具有同加性，同乘性，同幂性
线性同余方程
给定整数 \(a,b,m\)，求一个整数 \(x\) 满足：\(ax\equiv b \pmod m\)，或者判定无解。该方程等价于 \(ax-b\) 是 \(m\) 的整数倍，设为 \(-y\) 倍。
则：\(ax-b=m*(-y)\)。移向得：\(ax+my=b\)
这样就转化为了解线性不定方程，当 \(gcd(a,m)|b\) 时有解。
求出一组通解后，快速转化为最小整数解的方法：
设 \(ans=b/gcd(a,m)*x_0\)，\(s=m/gcd(a,m)\)，则最小整数解为 \(ans=(ans\%s+s)%s\)。（仅用于同余方程，线性不定方程不适用）
ll exgcd(ll a, ll b, ll &x, ll &y){
if(b == 0){
x = 1, y = 0;
return a;
}
ll r = exgcd(b, a%b, x, y), t = x;
x = y;
y = t-(a/b)*y;
return r;
}
int cal(ll a, ll b, ll m, ll &x, ll &y){
ll g = exgcd(a, m, x, y);
if(b%g) return 0;
ll k = b/g, t = m/g;
x = (x*k%t+t)%t;
return 1;
}
欧拉定理 费马小定理
欧拉定理
欧拉定理：若正整数 \(a,n\) 互质，则 \(a^{\phi(n)}\equiv 1(mod\text{ }n)\)，其中 \(\phi(n)\) 为欧拉函数。
费马小定理
费马小定理：若正整数 \(a,n\) 互质，且 \(n\) 是质数，则 \(a^{n-1}\equiv 1(mod\text{ }n)\)，其中 \(\phi(n)\) 为欧拉函数，两边同乘 \(a\)，则有 \(a^n\equiv n(mod\text{ }n)\)。
可以看出，费马小定理是欧拉定理在 \(n\) 为质数情况下的特例。
扩展欧拉定理
欧拉定理的推论：若正整数 \(a,n\) 互质，则对于任意正整数 \(b\)，满足：\(a^b\equiv a^{b\text{ }mod\text{ }\phi(m)}(mod\text{ }m)\)
证明：设 \(b=q\times \phi(m)+r\)，则 \(a^r = a^{b\text{ }mod\text{ }\phi(m)}\)。因为 \(a^{\phi(m)} \equiv 1(mod\text{ }m)\)，所以 \(a^b\equiv a^{\phi(m)^q}\times a^r \equiv 1^q\times a^r \equiv a^r \equiv a^{\phi(m)}(mod\text{ }m)\)
再次扩展：若正整数 \(a,m\) 不一定互质，且对于任意正整数 \(b\) 且 \(b\ge \phi(m)\) 时，满足：\(a^b\equiv a^{b\text{ }mod\text{ }\phi(m)+phi(m)}(mod\text{ }m)\)
而对于 \(b < \phi(m)\)，则满足不能降幂，\(a^b \equiv a^b\)
void eular(){
a = read(), m = read();
p = phi(m);
while(1){//高精度b
c = getchar();
if(c < '0' || c > '9') break;
b = b*10+c-'0';
if(!f) tmp = tmp*10+c-'0';
if(!f && tmp > p) f = 1;
b %= p;
}
if(f) b += p;
if(!f) b = tmp;
write(qpow(a, b, m));
}
乘法逆元
逆元定义：模意义下的倒数。
若整数 \(b,m\) 互质，并且 \(b|a\)，则存在一个整数 \(x\)，使得 \(\frac{a}{b} \equiv a*x \pmod m\)，称 \(x\) 为 \(b\) 的模 \(m\) 乘法逆元，记为 \(b^{-1}\pmod m\)。
因为 \(\frac{a}{b}\equiv a/b*b*b^{-1}\equiv a*b^{-1} \pmod m\)，所以 \(b*b^{-1} \equiv 1 \pmod m\)，该过程需要满足 \(b|a\)。求 \(b^{-1}\) 相当于求 \(bx\equiv 1\pmod m\) 的解。该同余方程需要满足 \(b,m\) 互质。
所以 \(b|a\) 和 \(b,m\) 互质是求乘法逆元的两个前提。
特别的，如果 \(m\) 是质数，可通过费马小定理得出 \(b\) 的乘法逆元即为 \(b^{m-2}\)，通过快速幂直接求即可。
ll inv(ll b, ll p, ll op){
if(!op){ // p为质数
return qpow(b, p-2);
}
else{ //正常情况(b和p互质)
exgcd(b, p, x, y);
x = (x%p+p)%p; //因为是互质的b与p求最大公约数，所以exgcd返回值为1，原解线性同余方程所需乘的系数就不需要了
return x;
}
}
线性求 \(1\sim p-1\) 的模 \(p\) 下的逆元(\(p\) 为质数)：
\(inv[1] = 1,p=k*i+r(k=p/i,r=p\%i)\)
\(k*i\equiv -r\pmod p\)
同时乘 \(i^{-1}r^{-1}\)
\(kr^{-1}\equiv -i^{-1}\pmod p\)
\(i^{-1} \equiv -k(p\%i)^{-1}\pmod p\)
因为 \(0\) 和 \(p\) 在模 \(p\) 意义下可互换
\(i^{-1}\equiv (p-k)(p\%i)^{-1}\pmod p\)
\(i^{-1}\equiv (p-p/i)(p\%i)^{-1}\pmod p\)
所以 \(inv[i]=(p-p/i)inv[p\%i]\%p\)
void allinv(){
inv[1] = 1;
for(ll i = 2;i <= n;i++) inv[i] = (p-p/i)*inv[p%i]%p;
}
中国剩余定理及其扩展
中国剩余定理
中国剩余定理（简称 crt），主要用于解决 \(\left\{\begin{matrix}x \equiv a_1 (mod\text{ }m_1)\\x \equiv a_2 (mod\text{ }m_2)\\...\\...\\x \equiv a_n (mod\text{ }m_n)\end{matrix}\right.\) 且 \(m_1,m_2,...,m_n\) 互质类的问题。
设 \(M=\Pi _{i=1}^nm_i\)，\(M_i=M/m_i\)，\(M_i^{-1}\times M_i\equiv 1\pmod {m_i}\)，即 \(M_i^{-1}\) 为 \(M_i\) 的乘法逆元，模意义下的倒数。因为 \(\left\{\begin{matrix}a_jM_jM_j^{-1}\%m_i=0,j\ne i\\a_jM_jM_j^{-1}\%m_i=a_i,j= i\end{matrix}\right.\) ，所以可以求出一个解为 \(x = \Sigma_{i=1}^na_iM_iM_i^{-1}\)。由此，任意解 \(x_0即为x+n\times M\)，最小正整数解\(x_{min}\)即为\(x\%M\)。由于需要求逆元，所以需要 \(exgcd\) 算法辅助。
void crt(){
for(ll i = 1;i <= n;i++){
exgcd(sum/m[i], m[i], x, y);
ans = (ans+a[i]*x*sum/m[i])%sum; //sum/m[i]即为M[i]，x即为逆元
}
write((ans+sum)%sum); //防止负数
}
扩展中国剩余定理
当模数\(m_1,m_2,...,m_n\)不互质时，一般的crt无法求解，就需要用到扩展中国剩余定理（excrt）了。
数学推导
举个例子：
\(\left\{\begin{matrix} x\equiv2\pmod 4 \\x\equiv4\pmod 6 \end{matrix}\right.\) \(\Rightarrow\) \(x\equiv10\pmod {12}\)
\(\left\{\begin{matrix} x\equiv4\pmod 6 \\x\equiv3\pmod 5 \end{matrix}\right.\) \(\Rightarrow\) \(x\equiv 28\pmod {30}\)
\(\left\{\begin{matrix} x\equiv2\pmod 4 \\x\equiv3\pmod 6 \end{matrix}\right.\) \(\Rightarrow\) \(\phi\)
可以粗略地总结出几个规律：
新方程与原方程具有同样的形式；
新方程的模数，是之前两个模数的lcm；
方程可能存在无解的情况
形式化的，对于
\(\left\{\begin{matrix} x\equiv {r_1}\pmod {m_1} \\x\equiv {r_2}\pmod {m_2} \end{matrix}\right.\) \(\Rightarrow\) \(a=k_1m_1+r_1=k_2m_2+r_2\) \(\Rightarrow\) \(k_1m_1-k_2m_2=r_2-r_1\)
其中\(m_1,m_2,r_1,r_2\)都为已知量，是一个典型的线性不定方程，可以直接用exgcd求解
对于如何求出一组\(k_1,k_2\)，设\(g=\gcd(m_1,m_2),p_1=\frac{m_1}{g},p_2=\frac{m_2}{g}\)，显然\(p_1,p_2\)互质
上面的式子可改写成\(k_1p_1-k_2p_2=\frac{r_2-r_1}{g}\)，因为仅当\(g|(r_2-r_1)\)时才有解，所以右边为整数，又\(p_1,p_2\)互质，求解这个方程显然能用exgcd
设解出的两个解为\(x_1,x_2\)，\(x_1p_1+x_2p_2=1\)，则有
\(\left\{\begin{matrix} k_1=\frac{r_2-r_1}{g}x_1 \\k_2=-\frac{r_2-r_1}{g}x_2 \end{matrix}\right.\)
于是一个特解\(x_0=r_1+k_1m_1=r_1+\frac{r_2-r_1}{g}x_1m_1\)
通解\(x=x_0+k\times lcm(m_1,m_2)\)。
程序实现
遍历每一个方程，若无解，直接输出，否则合并成一个方程，与下一个方程继续进行操作。
void excrt(){
ll M = a[1], x, y;
ans = b[1];
for(ll i = 2;i <= n;i++){
ll g = exgcd(M, a[i], x, y), c = ((b[i]-ans)%a[i]+a[i])%a[i];
if(c%g){
puts("-1");
return;
}
x = c/g*x%a[i];
ans += x*M;
M = a[i]/g*M;
ans = (ans%M+M)%M;
}
}
整除分块与简单数列
整除分块
推论 1：对于 \(i\in [x,k/(k/x)]\)，\(k/i\) 的值都相等。
如 \(k=15,x=9\)，则区间 \([9,15/(15/9)]\)，即对于 \(\forall i\in[9,15]\)，\(k/i\) 的值都相等。
推论 2：对于 \(i\in[1,k]\)，\(k/i\) 至多只有 \(2\sqrt k\) 个不同的值。
如当 \(k=15\) 时，\(k/i\) 的结果至多只有 \(6\) 个。
可以通过推论 1 把 \(1-k\) 分成若干个块，加速遍历；通过推论 2 可以推算时间复杂度。
简单数列
整除分块经常配合等差、等比数列求解，这里给出几个公式：
等差数列
对于数列 \(a_1,a_2,....,a_n\)，设公差为 \(d\)，设 \(S_n\) 表示等差数列的前 \(n\) 项和。
\(S_n=\frac{1}{2}n(a_1+a_n)=\frac{d}{2}n^2+(a_1-\frac{d}{2})n\)
等比数列
对于数列 \(a_1,a_2,....,a_n\)，设公比为 \(q\)，设 \(S_n\) 表示等比数列的前 \(n\) 项和。
\(S_n=\frac{a_1-a_nq}{1-q}=\frac{a_1(1-)}{1-q}(q\ne1)\)
BSGS 算法及其扩展
BSGS
BSGS，即 baby step giant step（大步小步算法），是一种通过分块思想求解高次同余方程的算法。
普通BSGS只能求解当\(a,p\)互质时的问题。
因为\(a,p\)互质，所以\(a^x\equiv a^{x\bmod \phi(p)}\pmod p\)，因为\(\phi(p)<p\)所以\(x<p\)，所以可以设\(x=i*t-j\)，其中\(t=\lceil{\sqrt p}\rceil\)，\(0 \le j \le t-1\)，则方程变为\(a^{i*t-j}\equiv b \pmod p\)，即\((a^t)^i\equiv b*a^j \pmod p\)。
对于所有的\(j \in [0,t-1]\)，把\(b*a^j \bmod p\)插入一个哈希表。
枚举\(i\)的所有可能取值，即\(i \in [0,t]\)，计算出\((a^t)^i \bmod p\)，在哈希表中查找是否存在对应的\(j\)，更新答案即可。
时间复杂度为\(O(\sqrt p)\)。
ll BSGS(ll x, ll y, ll p){
map <ll, ll> mp;
mp.clear();
y %= p;
ll t = ceil(1.0*sqrt(p));
for(ll j = 0;j < t;j++) mp[y*qpow(x, j, p)%p] = j;
x = qpow(x, t, p);
if(x == 0){
if(y == 0) return 1;
return -1;
}
for(ll i = 0;i <= t;i++){
ll val = qpow(x, i, p);
if(mp.find(val) != mp.end() && i*t-mp[val] >= 0) return i*t-mp[val];
}
return -1;
}
扩展BSGS
狄利克雷卷积&莫比乌斯反演
组合数学
排列和组合
排列
从 \(n\) 个不同元素中取出 \(m\) 个元素，按照一定的顺序排成一列，叫做从 \(n\) 个元素中取出 \(m\) 个元素的一个排列。相同的排列元素和排列顺序都相同。
从 \(n\) 个不同元素中取出 \(m(m \le n)\) 个元素的所有排列个数，叫做排列数，记为 \(A_n^m\)。
\(A_n^m=n(n-1)(n-2)...(n-m+1)=\frac{n!}{(n-m)!}\)。
特别的，当 \(m=n\) 时，叫做全排列。\(A_n^n=n!\)。规定 \(0!=1\)。
组合
从 \(n\) 个不同元素中取出 \(m\) 个元素并成一组，不讲究顺序，叫做从 \(n\) 个元素中取出 \(m\) 个元素的一个组合。相同的排列元素和排列顺序都相同。
从 \(n\) 个不同元素中取出 \(m(m \le n)\) 个元素的所有组合个数，叫做组合数，记为 \(C_n^m\)。
\(C_n^m=\frac{n(n-1)(n-2)...(n-m+1)}{m!}=\frac{n!}{(n-m)!\ \cdot \ m!}\)。
性质：
1.\(C_n^m=C_n^{n-m}\)（规定 \(C_n^0=1\)）
2.\(C_n^m=C_{n-1}^m+C_{n-1}^{m-1}\)（相当于对于 \(n\) 号元素取或不取两种情况）
3.\(C_n^0+C_n^1+...+c_n^n=2^n\)（即为每个元素都选不不选两种选择，故 \(2^n\)）
组合数的求法
根据上式 2，可通过递归求出 \(0\le y\le x\le n\) 的所有组合数 \(C_x^y\), 时间复杂度 \(O(n^2)\)。
void solve1(){
for(int i = 0;i <= x;i++) c[i][0] = c[i][i] = 1;
for(int i = 1;i <= x;i++){
for(int j = 1;j <= i/2;j++){
c[i][j] = c[i][i-j] = c[i-1][j-1]+c[i-1][j];
}
}
}
若要求出对 \(C_n^m\) 对 \(p\) 取模后的值，并且 \(1-n\) 都存在乘法逆元，则可以先计算分子 \(n!\bmod p\), 再计算分母 \(m!(n-m)!\bmod p\) 的逆元，乘起来后得到 \(C_n^m\bmod p\), 复杂度为 \(O(n)\)。
若在计算阶乘过程中，把 \(0\le k\le n\) 的每个 \(k!\bmod p\) 及其逆元分别保存在数组中，则可以通过 \(O(n\log n)\) 的预处理后 \(O(1)\) 获得 \(C_x^y \bmod p=jc[x]*jc_inv[y]*jc_inv[x-y]\bmod p\)。
若需要对 \(C_n^m\) 进行高精度计算，为避免除法，可将分子分母的质因数的指数保存下来，相减（把分母除去）后把剩余质因数乘起来，时间复杂度 \(O(n\log n)\)。
二项式定理
\((a+b)^n=\Sigma_{i=0}^nC_n^ia^{n-i}b^i\)
二项式定理与杨辉三角
杨辉三角：
1
1 1
1 2 1
1 3 3 1
1 4 6 4 1
1 5 10 10 5 1
....
\((a+b)^0 = 1\)
\((a+b)^1 = a+b\)
\((a+b)^2=a^2+2ab+b^2\)
....
\((a+b)^5=a^5+5a^4b+10a^3b^2+10a^2b^3+5ab^4+b^5\)
多重集
Lucas 定理及其扩展
Lucas 定理
若 \(p\) 是质数，则对于任意整数 \(1\le m\le n\)，有 \(C_n^m \equiv C_{n\bmod p}^{m \bmod p}*C_{n/p}^{m/p}\pmod p\)。
可以理解为把 \(n,m\) 表示为 \(p\) 进制数，对 \(p\) 进制下的每一位分别计算组合数，最后再乘起来。
Lucas 定理主要用于求解 \(n,m\) 很大，\(p\) 较小的 \(C_n^m\bmod p\) 问题。
经过递归求解，时间复杂度为 \(O(p+\log_pn*log_2p)\)。
ll c(ll n, ll m){
if(m > n) return 0;
return pr[n]*qpow(pr[n-m], p-2)%p*qpow(pr[m], p-2)%p; //pr为阶乘
}
ll lucas(ll n, ll m){
if(m == 0) return 1;
return c(n%p, m%p)*lucas(n/p, m/p)%p;
}
若\(p\)较大时，预处理阶乘复杂度\(O(p)\)会TLE。此时可以改变求组合数\(C\)的方式，通过\(C_n^m=\frac{n*(n-1)*...*(n-m+1)}{m!}\)，逐项求出分子和分母模\(p\)的值，然后用快速幂求乘法逆元就可以了。
扩展 Lucas 定理
对于 \(p\) 是和数的情况，就需要用到扩展 Lucas 定理。
将 \(p\) 分解为质因数的形式，即 \(p=p_1^{k_1}*p_2^{k_2}*...*p_n^{k_n}\)。列出同余方程式：
\(\left\{\begin{matrix}ans \equiv c_1 \pmod {{p_1}^{k_1}}\\ans \equiv c_2 \pmod {{p_2}^{k_2}} \\...\\ans \equiv c_q \pmod {{p_q}^{k_q}}\end{matrix}\right.\)
其中\(c_1...c_q\)为每个\(C_n^m\%p_i^{k_i}\)求出的答案
1.因为\(p_i^{k_i}\)为质因数分解后的结果，所以满足互质，可以直接用中国剩余定理合并
2.通过求出\(n!\%p_i^{k_i},m!\%p_i^{k_i},(n-m)!\%p_i^{k_i}\)再通过逆元求出\(C_n^m\%p_i^{k_i}\)
3.如何求出\(n!\%p_i^{k_i}\)：
以\(n=19,p_i=3,k_i=2\)举例：
\(n!=1*2*3*4*5*6*7*8*9*10*11*12*13*14*15*16*17*18*19\\ \ \ \ =(1*2*4*5*7*8*10*11*13*14*16*17*19)*3^6*(1*2*3*4*5*6)\)
根据这个例子发现，求解\(n!\)可以分为三部分：
第一部分是\(p_i\)的幂部分\(p_i^{\lfloor{\frac{n}{p_i}}\rfloor}\)，即\(3^6\)。
第二部分是一个新的新的阶乘。
第三部分是剩下的数。发现这部分在\(\bmod p_i^{k_i}\)下是以\(p_i^{k_i}\)为周期的，即
\((1*2*4*5*7*8)\equiv (10*11*13*14*16*17)\pmod {p_i^{k_i}}\)，所以只需求\(p_i^{k_i}\)的长度即可。还剩下孤立的数长度不会超过\(p_i^{k_i}\)，暴力求解即可。
4.求出的\(m!\%p_i^{k_i}\)和\((n-m)!\%p_i^{k_i}\)有可能与\(p_i^{k_i}\)不互质，无法求逆元
所以要将\(m!\%p_i^{k_i}\)和\((n-m)!\%p_i^{k_i}\)中的质因子\(p_i\)全部除去。求出逆元后再全部乘回去
计算\(n!\)中质因子\(p\)的公式为\(x=\lfloor{\frac{n}{p}}\rfloor+\lfloor{\frac{n}{p^2}}\rfloor+\lfloor{\frac{n}{p^3}}\rfloor+...\)，递归式为\(f(n)=f(\lfloor{\frac{n}{p}}\rfloor)+\lfloor{\frac{n}{p}}\rfloor\)。
至此，整个问题便可以完美求解，时间复杂度为\(O(p+\log n\log p)\)。
ll cal(ll n, ll pi, ll pk){
if(n == 0) return 1;
ll ans = 1;
if(n/pk){
for(ll i = 1;i <= pk;i++) if(i%pi) ans = ans*i%pk;
ans = qpow(ans, floor(n/pk), pk);
}
for(ll i = 1;i <= n%pk;i++) if(i%pi) ans = ans*i%pk;
return ans*cal(n/pi, pi, pk)%pk;
}
ll c(ll n, ll m, ll p, ll pi, ll pk){
if(m > n) return 0;
ll a = cal(n, pi, pk), b = cal(m, pi, pk), c = cal(n-m, pi, pk), k = 0, ans;
for(ll i = n;i;i /= pi) k += i/pi;
for(ll i = m;i;i /= pi) k -= i/pi;
for(ll i = n-m;i;i /= pi) k -= i/pi;
ans = a*inv(b, pk)%pk*inv(c, pk)%pk*qpow(pi, k, pk)%pk; //a,b,c已经在cal()中除过了
ans = ans*(p/pk)%p*inv(p/pk, pk)%p;
return ans;
}
ll lucas(ll n, ll m, ll p){
ll x = p, pk, ans = 0;
for(ll i = 2;i <= p;i++){
if(x%i) continue;
pk = 1;
while(x%i == 0) x /= i, pk *= i;
ans = (ans+c(n, m, p, i, pk))%p;
}
return ans;
}
Catalan 数
容斥原理
线性代数
矩阵快速幂
矩阵
矩阵在程序中通常以二维数组的形式出现，如\(m\)行\(n\)列的矩阵\(A\)：
\(\begin{bmatrix} A_{11} & A_{12} & ... & A_{1n}\\ A_{21} & A_{22} & ... & A_{2n}\\ ... & ... & ... & ... \\ A_{m1} & A_{m2} & ... & A_{mn} \end{bmatrix}\)
矩阵乘法
矩阵乘法要求矩阵\(A\)的列数必须和矩阵\(B\)的行数相等。即\(A_{mn}*B_{np}=C_{mp}\)（\(mn，np\)中间两位需一样，结果为两边的\(mp\)})。
对于矩阵\(C_{ij}\)，是\(A\)的第\(i\)行分别乘以\(B\)的第\(j\)列的和。计算三重循环即可
for(int i = 1;i <= m;i++)
for(int j = 1;j <= n;j++)
for(int k = 1;k <= p;k++)
c[i][j] += a[i][k]*b[k][j]
矩阵快速幂
以斐波那契数列为例：
\(f[n]=1*f[n-1]+1*f[n-2]\)
\(f[n-1]=1*f[n-1]+0*f[n-2]\)
可以转化为：
\(\begin{bmatrix} 1 & 1\\ 1 & 0 \end{bmatrix}\) \(\times\) \(\begin{bmatrix} f[n-1] \\ f[n-2] \end{bmatrix}\) = \(\begin{bmatrix} f[n] \\ f[n-1] \end{bmatrix}\)。
展开后为\(f[n]=f[i-1]+f[i-2],f[i-1]=f[i-1]\)。
使用矩阵快速幂求解的问题需满足\(K\)阶常系数线性递推关系。
以上面为例，\(\begin{bmatrix} f[n-1] \\ f[n-2] \end{bmatrix}\) 叫做状态矩阵，\(\begin{bmatrix} 1 & 1\\ 1 & 0 \end{bmatrix}\) 叫做转移矩阵，通常以\(A\)表示。
由上式可推得\(\begin{bmatrix} f[n] \\ f[n-1] \end{bmatrix}\) = \(A^{n-2}*\) \(\begin{bmatrix} 1 \\ 1 \end{bmatrix}\)。
所以求解\(f_n\)的关键就在于求矩阵\(A\)的\(n-2\)次方，用的方法就是快速幂。因为所求对象是矩阵，因此称为矩阵快速幂。矩阵用来加速递推，效率很高。对比普通快速幂的求法，用\(r\)来保存结果，并初始化。
void mul(){
for(ll i = 1;i <= n;i++){
for(ll j = 1;j <= p;j++){
t[i][j] = 0;
for(ll k = 1;k <= m;k++) t[i][j] = (t[i][j]+a[i][k]*r[k][j])%M;
}
}
for(ll i = 1;i <= n;i++) for(ll j = 1;j <= p;j++) r[i][j] = t[i][j];
}
void mulself(){
for(ll i = 1;i <= n;i++){
for(ll j = 1;j <= n;j++){
t[i][j] = 0;
for(ll k = 1;k <= m;k++) t[i][j] = (t[i][j]+a[i][k]*a[k][j])%M;
}
}
for(ll i = 1;i <= n;i++) for(ll j = 1;j <= n;j++) a[i][j] = t[i][j];
}
void qpow(ll x){
while(x){
if(x&1) mul();
mulself();
x >>= 1;
}
}
高斯消元
概念
高斯消元是一种求解线性方程组的方法。一个由\(M\)个\(N\)元一次方程共同构成的线性方程组，其中所有系数可以写成一个\(M\)行\(N\)列的矩阵，再加上方程等号右侧的常数，可以写成一个\(M\)行\(N+1\)列的增广矩阵。
如\(\left\{\begin{matrix} x_1+2x_2-x_3=-6 \\2x_1+x_2-3x_3=-9 \\-x_1-x_2+2x_3=7 \end{matrix}\right.\) $\Rightarrow $ \(\begin{bmatrix} 1 & 2 & -1 & -6\\ 2 & 1 & -3 & -9\\ -1 & -1 & 2 & 7 \end{bmatrix}\)。
解题思路即为把方程组列出，用高斯消元求解即可。
消元过程
求解操作：
1.用一个非\(0\)的数乘某一行
2.把其中一行的若干倍加到另一行上
3.交换两行的位置
上述操作被称为矩阵的“初等行变换”。
仍以上式举例：
\(\begin{bmatrix} 1 & 2 & -1 & -6\\ 2 & 1 & -3 & -9\\ -1 & -1 & 2 & 7 \end{bmatrix}\) $\Rightarrow $ \(\begin{bmatrix} 1 & 2 & -1 & -6\\ 0 & 1 & 1 & 1\\ 0 & 0 & 1 & 3 \end{bmatrix}\) $\Rightarrow $ \(\left\{\begin{matrix} x_1+2x_2-x_3=-6 \\x_2+x_3=1 \\x_3=3 \end{matrix}\right.\)
上图中间的矩阵被称为“阶梯型矩阵”。
\(\begin{bmatrix} 1 & 2 & -1 & -6\\ 0 & 1 & 1 & 1\\ 0 & 0 & 1 & 3 \end{bmatrix}\) \(\Rightarrow\) \(\begin{bmatrix} 1 & 0 & 0 & 1\\ 0 & 1 & 0 & -2\\ 0 & 0 & 1 & 3 \end{bmatrix}\) \(\Rightarrow\) \(\left\{\begin{matrix} x_1= 1 \\x_2 = -2 \\x_3 = 3 \end{matrix}\right.\)。
最后得到的矩阵被称为“简化阶梯型矩阵”，他的系数矩阵部分是一个“对角矩阵”。
高斯消元算法的思想就是对每个未知量\(x_i\)，找到一个\(x_i\)的系数非零，但\(x_1\sim x_{i-1}\)的系数都为零的方程，然后用初等行变换把其他方程的\(x_i\)的系数全部消成零。
特殊情况
上述的例子是一种比较理想的情况，事实上，求解过程中还会有很多特殊情况：
1.出现\(0=d\ (d \ne 0)\)的方程，此时方程组无解
2.找不到一个\(x_i\)的系数非零，但\(x_1\sim x_{i-1}\)的系数都为零的方程，例如：
\(\left\{\begin{matrix} x_1+2x_2-x_3=3 \\2x_1+4x_2-8x_3=0 \\-x_1-2x_2+6x_3=2 \end{matrix}\right.\) \(\Rightarrow\) \(\begin{bmatrix} 1 & 2 & 0 & 4\\ 0 & 0 & 1 & 1\\ 0 & 0 & 0 & 0 \end{bmatrix}\)。
此时找不到\(x_2\)的系数非零，但\(x_1\)的系数为零的方程，方程组的解可写作：\(\left\{\begin{matrix} x_1=4-2x_2 \\x_3=1 \end{matrix}\right.\)。
\(x_2\)取任意值，都可以计算出一个相应的\(x_1\)。此时方程组有无限多的解。我们把\(x_1,x_3\)称作主元，\(x_2\)称作自由元。
仔细分析可以发现，对于每个主元，简化阶梯形矩阵有且仅有一个位置\((i,j)\)，满足该主元的系数非零，其他都为零。
综上：
若存在系数全为零，常数不为零的行，则方程组无解
若系数不全为零的行恰好\(n\)个，则说明主元有\(n\)个，方程组有唯一解
若系数不全为零的行有\(k < n\)个，则说明主元有\(k\)个，自由元有\(n-k\)个，方程组有无限多的解。
注意事项
高斯消元法的细节很多，包括并不限于：
1.交换时不要只枚举\(i\)到\(n\)，应该从\(1\)到\(n\)全部检查一遍（因为可能存在0）
2.用eps判0
3.fabs()
4.交换，消元时不要忘了循环到\(n+1\)
5.除之前先判断除数是不是0
bool cmp(int i, int j){
if(fabs(a[j][i]) > fabs(a[i][i])) return 1; // a[j][i]为第j行的第i个，和a[i][i]在同一列
else if(fabs(a[j][i]) < fabs(a[i][i])) return 0;
for(int k = i+1;k <= n;k++){ // 尽可能让绝对值大的在上面
if(fabs(a[j][k]) > fabs(a[i][k])) return 1;
else if(fabs(a[j][k]) < fabs(a[i][k])) return 0;
}
return 0;
}
void gauss(){
for(int i = 1;i <= n;i++){
for(int j = 1;j <= n;j++){
if(fabs(a[j][j]) > eps && j <= i) continue; // 帮助判断无解或无限解
if(cmp(i, j)) for(int k = 1;k <= n+1;k++) swap(a[j][k], a[i][k]);
// 循环找第i列中绝对值最大的那一行
}
if(fabs(a[i][i]) < eps) continue;
for(int j = 1;j <= n;j++){
if(i == j || !a[j][i]) continue;
double r = a[i][i]/a[j][i]; // 上面的交换就是在尽可能让r变大，减小误差
for(int k = 1;k <= n+1;k++) a[j][k] = a[i][k]-a[j][k]*r; // 用a[i][i]消去所有行的a[j][i]
}
}
}
void solve(){
gauss();
for(int i = 1;i <= n;i++){
if(fabs(a[i][i]) < eps && fabs(a[i][n+1]) > eps){
puts("No Solution"); // 系数为零但常数非零
return;
}
}
for(int i = 1;i <= n;i++){
if(fabs(a[i][i]) < eps && fabs(a[i][n+1]) < eps){
puts("Countless Solution"); // 系数常数都为零
return;
}
}
for(int i = 1;i <= n;i++){
double ans = a[i][n+1]/a[i][i];
if(fabs(ans) < eps) ans = 0.00;
printf("x%d=%.2lf\n", i, ans);
}
}
线性空间
概率与期望
基本概念
在概率论中，把一个随机试验的某种可能结果称为“样本点”，把所有可能结果构成的集合称为“样本空间”。
在一个给定的样本空间中，随机事件就是样本空间的子集，即由若干个样本点构成的集合。
随机变量是试验结果样本点的实数化。
我们一般以\(P(A)\)表示随机事件A的发生几率，即随机事件A发生的可能性，是一个0~1的实数。
假设随机事件\(A\)属于样本空间\(\Omega\)，则
1.\(0\le P(A)\le 1\)
2.\(P(\Omega) = 1\)
3.假设若干随机事件\(A_1,A_2,A_3\)两两互斥，则\(\Sigma P(A_i)=P(\)
贝叶斯公式
0/1分数规划
0/1分数规划是指，给定整数\(a_1,a_2,..,a_n\)和\(b_1,b_2,...,b_n\)，求一组解\(x_i\ (1\le i\le n,\ x_i=0\ 或\ 1)\)。使得\(\frac{\Sigma_{i=1}^na_i*x_i}{\Sigma_{i=1}^nb_i*x_i}\)取得最大值。
简便来说，就是选出\(n\)对\(a_i,b_i\)，使选出的\(\frac{sum_a}{sum_b}\)最大。
我们猜测一个值\(L\)，然后考虑是否存在一组解，使\(\Sigma_{i=1}^n(a_i-L*b_i)*x_i\ge 0\)，变形后得\(\exists\left \{x_1,x_2,...,x_n \right \}，使得\frac{\Sigma_{i=1}^na_i*x_i}{\Sigma_{i=1}^nb_i*x_i}\ge L\)，也就是说，\(L\)比我们求的最大值要小。
同理，如果对于任意一组解，使\(\Sigma_{i=1}^n(a_i-L*b_i)*x_i<0\)，即\(\forall\left \{x_1,x_2,...,x_n \right \}，都有\frac{\Sigma_{i=1}^na_i*x_i}{\Sigma_{i=1}^nb_i*x_i}<L\)，也就是说，\(L\)比我们求的最大值要大。
上面的过程，与二分答案很相似，即可以通过猜测\(L\)得出答案应该更大或更小，也就是说\(L\)对于解的存在性具有单调性。
对于如何判定\(\Sigma_{i=1}^n(a_i-L*b_i)*x_i\ge 0\)，显然，只需算出所有\(a_i-L*b_i\)，将所有正数加起来就可以判断了。
综上，该问题可以通过实数二分解决。
bool check(double x){
double maxn = 0;
for(int i = 1;i <= n;i++){
t[i] = a[i]-1.0*x*b[i];
maxn = max(maxn, t[i]);
}
if(maxn > 0) return 1;
return 0;
}
void solve(){
double l = 0, r = maxr, mid;
while(r-l > eps){
mid = (l+r)/2;
if(check(mid)) l = mid;
else r = mid;
}
cout << l << endl;
}
博弈论
基础概念
公平组合游戏ICG
若一个游戏满足
1.由两名玩家交替操作
2.在游戏进程的任意时刻，可以执行的合法行动与轮到哪位玩家无关
3.不能行动的玩家判负
则被成为一个公平组合游戏
有向图游戏
给定一个有向无环图（DAG），图中有唯一的起点，起点上有一个棋子，两名玩家交替将这枚棋子沿有向边进行移动，一次一步，无法移动者判负，该游戏被称为有向图游戏。
任何一个公平组合游戏都可以转化为有向图游戏，有些题目中可以以图为突破点。
Mex运算
设\(mex(S)\)为求出不属于集合\(S\)的最小非负整数运算。即\(mex(S)=min\left \{ x \right \}\ (x \notin S, x \in N)\)。
如\(mex(1,2,3)=0,mex(0,1,3)=2,mex()=0\)。
Nim博弈（尼姆博弈）
这里介绍一个基础的博弈游戏，后面会以此为例。
定义：\(n\)堆物品，每堆物品的个数任意，两人轮流取，每次取某堆中不少于1个，最后取完者胜。
结论：将每堆物品的数量全部异或起来，若值为0，则先手必败，否则先手必胜。
SG函数
在有向图游戏中，对于每个节点\(x\)，设从\(x\)出发有\(k\)条单向边，分别到达\(y_1,y_2,...,y_k\)，定义\(SG(x)\)为\(x\)的后继节点\(y_1,y_2,...y_k\)的\(SG\)函数值构成的集合再进行mex运算后的结果，即\(SG(x)=mex({SG(y_1),SG(y_2),...,SG(y_k)})\)
特别的，整个有向图游戏\(G\)的SG函数值被定义为有向图游戏起点\(s\)的SG函数值，即\(SG(G)=SG(s)\)。
例：以SG函数求解Nim博弈
我们可以将一个有\(x\)个物品的堆视为节点\(x\)，则当且仅当\(y<x\)时，节点\(x\)可以到达\(y\)。
这样即可把由\(n\)个堆组成的Nim游戏转化成\(n\)个有向图游戏，再用SG函数求解即可。
定理
有向图游戏的某个局面必胜，当且仅当该局面对应节点的SG函数值大于0。
有向图游戏的某个局面必败，当且仅当该局面对应节点的SG函数值等于0。
有向图的和
设\(G_1,G_2,...,G_m\)是\(m\)个有向图游戏。定义有向图游戏\(G\)，其行动规则是任选某个有向图游戏\(G_i\)，并在其上行动一步。\(G\)被称为有向图游戏\(G_1,G_2,...,G_m\)的和。有向图游戏的和的SG函数值等于它包含的各个子游戏SG函数值的异或和，即:
\(SG(G)=SG(G_1)\ xor\ SG(G_2)\ xor\ ...\ xor\ ...\ SG(G_m)\)。
应用
1.可选步数为1~m的连续整数，直接取模即可，\(SG(x) = x\%(m+1)\);
2.可选步数为任意步，\(SG(x) = x\);
3.可选步数为一系列不连续的数，用getsg()计算。
博弈论问题的一般解法：通过一些必胜必败情况推出上一步，通过问题情况推出公式后结合SG函数进行解答。
void getsg(){
for(int i = 1;i <= maxn;i++){
mst(mk, 0);
for(int j = 1;j <= k;j++){
mk[sg[i-b[j]]] = 1;
}
for(int j = 0;j <= maxn;j++){
if(!mk[j]){
sg[i] = j;
break;
}
}
}
}
参考资料
李煜东《算法竞赛进阶指南》
zzszsx的ppt
oi-wiki
众多cnblogs和csdn大神的博客
浙公网安备 33010602011771号