算法学习15：数论（高斯消元，组合数，卡特兰数）
文章目录
前言
提示：以下是本篇文章正文内容：
一、高斯消元
1.输入一个包含n个方程，n个未知数的线性方程组，求解这个方程组。
初等行列变换（线性代数）（最后得到一个行阶梯）
（原理）
1.把某一行乘以一个非零的数
2.交换某2行
3.把某行的若干倍加到另一行上去
（步骤实现） 枚举每一列：
1.找到每一列绝对值最大的一行
2.将该行移到最上面（性质2）
3.将该行的第一个数变为1（性质1）
4.将下面所有行的第c列都变成0（性质3）
输入一个包含n个方程，n个未知数的线性方程组。
求解这个方程组。
方程的解有三种可能：
1.有唯一解：输出n行，对应n个未知数的解
2.有无穷多解：infinite group solutions
3.无解：no sulution
#include <iostream>
#include <algorithm>
#include <cmath>
using namespace std;
const int N = 110;
const int eps = 1e-6;// 与0做判断用
// 0在计算机中存储可能是0.000000001（不用数有几个0，随便输入的）
// < eps表示为0， >eps表示不为0.
int n;
double a[N][N];// 系数矩阵 + 常数项
// 高斯消元
int gauss()
{
int c, r;// 列，行
for(c = 0, r = 0; c < n; c ++)
{
// 找到每一列绝对值最大的一行
// fabs：浮点数的绝对值
int t = r;
for(int i = r; i < n; i ++)
if(fabs(a[i][c]) > fabs(a[t][c]))
t = i;
// fabs(a[t][c]) == 0
// 这一列所有元素都为0，直接跳过
if(fabs(a[t][c]) < eps) continue;
// 将该行移到最上面（性质2）
// 将该行的第一个数变为1（性质1）
// 将下面所有行的第c列都变成0（性质3）
for(int i = c; i <= n; i ++) swap(a[t][i], a[r][i]);
for(int i = n; i >= c; i --) a[r][i] /= a[r][c];
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 ++)
// 列要多1，常数项矩阵
for(int j = 0; j < n + 1; j ++)
cin >> a[i][j];
int t = gauss();
//0：唯一解， 1：无穷多解， 2：无解
if(t == 0)
{
// 第n+1类存储的就是未知数的解
for(int i = 0; i < n; i ++)
printf("%.2lf\n", a[i][n]);
}
else if(t == 1) puts("Infinite group solutions");
else puts("No solution");
return 0;
}
二、组合数：
（1）给定n组询问，每组给定两个整数a，b，求Cab（从a个样品取出b个样品）mod（10^9 + 7）的值
求组合数1：1<n<10000，1<b<=a<2000（递推）
#include <iostream>
#include <algorithm>
using namespace std;
const int N = 2000 + 10, 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;// j = 0时，值为1
// 公式
else c[i][j] = (c[i - 1][j] + c[i - 1][j - 1]) % mod;
}
int main()
{
init();
int n;
scanf("%d", &n);
while(n --)
{
int a, b;
scanf("%d%d", &a, &b);
printf("%d\n", c[a][b]);
}
return 0;
}
组合数2：1<n<10000 1<b<=a<10^5（递推）
#include <iostream>
#include <algorithm>
using namespace std;
typedef long long LL;
const int N = 1e5 + 10, mod = 1e9 + 7;
int fact[N], infact[N];
// 计算：a^k % p
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()
{
fact[0] = infact[0] = 1;// 0! == 1
for(int i = 1; i < N; i ++)
{
fact[i] = (LL)fact[i - 1] * i % mod;
// qmi(i, mod - 2, mod)：i % mod 的逆元
infact[i] = (LL)infact[i - 1] * qmi(i, mod - 2, mod) % mod;
}
int n;
scanf("%d", &n);
while(n --)
{
int a, b;
scanf("%d%d", &a, &b);
// 除法 转换为逆元
printf("%d\n", (LL)fact[a] * infact[b] % mod * infact[a - b] % mod);
}
return 0;
}
（3）1<n<20 1<b<=a<10^18 1<p<10^5， 特点：a，b的范围特别大，卢卡斯定理（lucas）
#include <iostream>
#include <algorithm>
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;
a = (LL)a * a % p;
k >>= 1;
}
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;
}
// 卢卡斯定理
// 注意1：a，b的输入值可能比较大
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);
}
int main()
{
int n;
cin >> n;
while(n --)
{
LL a, b;
cin >> a >> b >> p;
cout << lucas(a, b) << endl;
}
return 0;
}
（4）1<b<=a<5000， 特点：只用算一组，但是结果特别大，要使用高精度计算
// 开启O2优化
// #pragma GCC optimize(2)
#include <iostream>
#include <algorithm>
#include <vector>
using namespace std;
typedef long long LL;
const int N = 5000 + 10;
int primes[N], cnt;// primes：存储1-a的质数 ， 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;
}
}
}
// 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;
// 1.筛素数
// 把1 - 5000 内的素数筛出来
get_primes(a);
// 质因数：24 = 2 * 2 * 2 * 3
// 2.求每一个质数出现的次数
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]);
for(int i = res.size() - 1; i >= 0; i --) printf("%d", res[i]);
puts("");
return 0;
}
三、卡特兰数：求有几种方案。
#include <iostream>
#include <algorithm>
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;
// 使用求组合数“最原始的公式”
for(int i = a; i > a - b; i --) res = (LL)res * i % mod;
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;
}
总结
提示：这里对文章进行总结：
💕💕💕
3355
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
