(资料图片仅供参考)
高斯消元
- 求解含有n个未知数,n个方程的多元线性方程组 O(n^3)
- 初等行变换:
- 某行乘以一个非零数
- 交换两行
- 某行加上另一行的若干倍
- 利用初等行变换将方程组化为上三角矩阵
- 解的情况:
- 完美阶梯型:唯一解
- 非完美阶梯型:
- 0 == 非0:无解
- 0 == 0:无穷解
- 步骤:
- 枚举每一列
- 找到这一列系数的绝对值最大的一行
- 将这一行与第一行交换
- 将改行的第一个数变成一(方程两边同乘某数)
- 把下面所有行的当前列的系数消成0(某行加上第一行的若干倍)
- 枚举每一列
- 代码:
const int N = 110;const double esp = 1e-6; //x < esp,则x = 0,否则 x != 0,由于浮点数精度问题,0可能是0.000001double a[N][N];int n;int gauss(){int r, c;for(r = 1, c = 1; c <= n; ++ c){//遍历每一列int t = r;for(int i = r; i <= n; ++ i) //找当前列的系数最大的行if(fabs(a[t][c]) < fabs(a[i][c]))t = i; //记录最大行if(fabs(a[t][c]) < esp) continue; //最大行系数为0说明该列系数均为0for(int i = c; i <= n + 1; ++ i) swap(a[r][i], a[t][i]);//否则交换第一行与系数最大行for(int i = n + 1; i >= c; -- i) a[r][i] /= a[r][c]; //将第一行当前列系数变为1for(int i = r + 1; i <= n; ++ i) //将下面所有行的当前列的系数变为0if(fabs(a[i][c]) > esp) //当前列系数非0for(int j = n + 1; 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 + 1]) > esp) return 0; //无解return 1; //无数解}for(int i = n; i >= 1; -- i) //从后往前求解for(int j = i + 1; j <= n; ++ j)a[i][n + 1] -= a[i][j] * a[j][n + 1];return 2; //唯一解}int main(){cin >> n;for(int i = 1; i <= n; ++ i)for(int j = 1; j <= n + 1; ++ j)cin >> a[i][j];int t = gauss();if(t == 1) cout << "Infinite group solutions";else if(t == 0) cout << "No solution";else if(t == 2){for(int i = 1; i <= n; ++ i) printf("%.2f\n", a[i][n + 1]);}return 0;}
求组合数
- \(C_a^b = \frac{a(a-1)(a-2)···(a-b+1)}{b!}=\frac{a!}{b!(a-b)!}\)
- 一般直接用公式求组合数容易超时
- 故可以通过预处理的方式,用递推式求出所有可能的组合数
- 递推式:\(C_a^b = C_{a-1}^{b-1} + C_{a-1}^{b}\) \(O(N^2)\)
- 当a,b <= 2e3,询问n = 1e4次时
- 代码:
//询问10000次//a,b <= 2000const int N = 2e3 + 10;const int mod = 1e9 + 7;int c[N][N];int 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 - 1] + c[i - 1][j]) % mod;}
- 当a,b <= 1e5,询问n = 1e4次,模数p = 1e9 + 7时
- 开二维数组会爆内存,于是利用公式,预处理出阶乘以及阶乘的逆\(O(Nlog(N))\)
- 代码:
typedef long long LL;const int N = 1e5 + 10;const int mod = 1e9 + 7;int fact[N], infact[N];int n;//快速幂int qmi(int a, int b, int p){ int res = 1; while(b){ if(b & 1) res = (LL)res * a % p; a = (LL)a * a % p; b >>= 1; } return res;}//预处理各阶乘以及阶乘的逆void init(){ fact[0] = infact[0] = 1; for(int i = 1; i < N; ++ i){ fact[i] = (LL)fact[i - 1] * i % mod; //递推式n! = (n - 1)! * n; infact[i] = (LL)infact[i - 1] * qmi(i, mod - 2, mod) % mod; //递推式(n!)^-1 = ((n - 1)!)^-1 * n^-1,取模意义下的逆就是逆元 }}int main(){ cin >> n; init(); while(n -- ){ int a, b; cin >> a >> b; cout << (LL)fact[a] * infact[b] % mod * infact[a - b] % mod << endl;//两个int最大值相乘不会爆longlong,但是三个相乘会爆 } return 0;}
- 当a,b <= 1e18,询问n = 2e5次,模数p不定时
- 两数相乘会爆longlong,利用卢卡斯定理预处理组合数 \(O(plog(N)log(p))\)
- 定理:\(C_a^b\equiv C_{a \mod p}^{b \mod p}·C_{a/p}^{b/p} \pmod p\)
- 代码:
typedef long long LL;int p;int n;//快速幂求逆元int qmi(int a, int b){ int res = 1; while(b){ if(b & 1) res = (LL)res * a % p; a = (LL)a * a % p; b >>= 1; } return res;}int C(int a, int b){ //求组合数C_a^b if(a < b) return 0; //可以不要 int x = 1, y = 1; //分子分母 for(int i = 0; i < b; ++ i){ x = (LL)x * (a - i) % p; //求分子 y = (LL)y * (i + 1) % p; //求分母 } return (LL)x * qmi(y, p - 2) % p;}//卢卡斯定理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(){ cin >> n; while(n -- ){ LL a, b; cin >> a >> b >> p; cout << lucas(a, b) << endl; } return 0;}
- 当a,b <= 5e5,询问一次,不取模时
- 此时组合数较大,需要用高精度存储结果
- 阶乘中素数的个数:
- \(v_p(n!) = \sum_{i=1}^{\infty} \left\lfloor \frac{n}{p^i} \right\rfloor\)
- 证明:将 \(n!\) 记为 \(1\times 2\times \cdots \times p\times \cdots \times 2p\times \cdots \times \lfloor n/p\rfloor p\times \cdots \times n\) 那么其中 \(p\) 的倍数有 \(p\times 2p\times \cdots \times \lfloor n/p\rfloor p=p^{\lfloor n/p\rfloor }\lfloor n/p\rfloor !\) 然后在 \(\lfloor n/p\rfloor !\) 中继续寻找 \(p\) 的倍数即可,这是一个递归的过程
- 代码:
const int N = 5e3 + 10;int primes[N], cnt;int st[N];int sum[N];//线性筛,1-n中质数就是n!中的质因子void get_primes(int a){ for(int i = 2; i <= a; ++ i){ if(!st[i]) primes[ ++ cnt] = i; for(int j = 1; primes[j] <= a / i; ++ j){ st[primes[j] * i] = 1; if(i % primes[j] == 0) break; } }}//求a!中质数p的个数int get(int a, int p){ int res = 0; while(a){ res += a / p; a /= p; } return res;}//高精度乘法vector
mul(vector &A, int b){ vector c; int t = 0; for(int i = 0; i < A.size() || t; ++ i){ if(i < A.size()) t += A[i] * b; c.push_back(t % 10); t /= 10; } while(c.size() > 1 && c.back() == 0) c.pop_back(); return c;}int main(){int a, b; cin >> a >> b; get_primes(a); for(int i = 1; i <= cnt; ++ i){ int p = primes[i]; sum[i] = get(a, p) - get(b, p) - get(a - b, p); //a!中质因子个数-b!中质因子个数-(a-b)!中质因子个数 } vector res; res.push_back(1); for(int i = 1; i <= cnt; ++ i) //累乘 for(int j = 1; j <= sum[i]; ++ j) res = mul(res, primes[i]); for(int i = res.size() - 1; i >= 0; -- i)//输出 cout << res[i]; return 0;}
卡特兰数
- \(C_{2n}^n-C_{2n}^{n-1}=\frac{1}{n+1}·C_{2n}^n\)
- 用于解决序列排列等问题
- AcWing 889. 满足条件的01序列 - AcWing
- 代码:
#include
using namespace std;typedef long long LL;const int N = 1e5 + 10;const int p = 1e9 + 7;int n;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;}int main(){ cin >> n; int res = 1; int a = 2 * n, b = n; int x = 1, y = 1; for(int i = 0; i < b; ++ i){ x = (LL)x * (a - i) % p; y = (LL)y * (i + 1) % p; } res = (LL)x * qmi(y, p - 2) % p; res = (LL)res * qmi(n + 1, p - 2) % p; cout << res; return 0;}