\(\color{black}高斯消元\)
题意
输入一个包含 \(n\) 个方程 \(n\) 个未知数的线性方程组。
方程组中的系数为实数。
求解这个方程组。
下图为一个包含 \(m\) 个方程 \(n\) 个未知数的线性方程组示例:
\[\left\{\begin{matrix}a_{1,1} \cdot x_1 + a_{1,2} \cdot x_2 + \dots + a_{1, n} \cdot x_n = b_1
\\a_{2,1} \cdot x_1 + a_{2,2} \cdot x_2 + \dots + a_{2, n} \cdot x_n = b_2
\\ \dots \dots \dots \dots \dots \dots \dots \dots \dots \dots \dots \dots\\a_{m,1} \cdot x_1 + a_{m,2} \cdot x_2 + \dots + a_{m, n} \cdot x_n = b_m
\end{matrix}\right.
\]
输入格式
第一行包含整数 \(n\)。
接下来 \(n\) 行,每行包含 \(n+1\) 个实数,表示一个方程的 \(n\) 个系数以及等号右侧的常数。
输出格式
如果给定线性方程组存在唯一解,则输出共 \(n\) 行,其中第 \(i\) 行输出第 \(i\) 个未知数的解,结果保留两位小数。
如果给定线性方程组存在无数解,则输出 Infinite group solutions。
如果给定线性方程组无解,则输出 No solution。
样例输入输出
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
转化
将方程组化成 1 个 \(n \cdot (n+1)\) 的系数矩阵,可以将其执行以下 3 种初等行列变换操作。
- 把某一行乘一个非零的数;
- 交换某 2 行;
- 把某行的若干倍加到另一行上去;
操作后,可以将其消成类似上三角的形式,即:
\[\left\{\begin{matrix}a_{1, 1} \cdot x_1 + a_{1, 2} \cdot x_2 + \dots + a_{1, n} \cdot x_n = b_1
\\a_{2,2} \cdot x_1 + a_{2,3} \cdot x_2 + \dots + a_{2, n} \cdot x_n = b_2
\\ \dots \dots \dots \dots \dots \dots \dots \dots \dots \dots
\\a_{n-1, n-1} \cdot x_{n-1} + a_{n-1, n} \cdot x_n = b_{n-1} \\x_n = b_n\end{matrix}\right.
\]
无解/有无穷多种解/有唯一解
- 有唯一解 : 消后是一个完美的阶梯型,即第 \(1\) 行有 \(n\) 个未知数,第 \(2\) 行有 \(n-1\) 个未知数,\(\cdots\cdots\) ,第 \(n\) 行有 \(1\) 个未知数。
- 无解 : \(0 = 非0\)
- 有无穷多组解 : \(0 = 0\)
求解
- 枚举每一列 \(c\),找到本列绝对值最大的一行.
- 将该行换到最上面
- 将该行的第 \(1\) 个数变成 \(1\),即将等式两边同时除以这个数
- 把第 \(c\) 列除了第 \(1\) 个数全部消成 \(0\)。
- 固定第 \(c\) 列,以后的操作将不影响第 \(1\) 到第 \(c\) 列的数值。
- 迭代操作结束后,将所有第 \(i\) 列第 \(i\) 行的数值消成 \(1\),再依次求解。
代码
int gauss(){
int c, r; // c 是列,r 是行
for(c = r = 0; c<n; c++){ // 枚举每一列
// 找到 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;
}
// 交换这两行
for(int i=c; i<=n; i++){
swap(a[t][i], a[r][i]);
}
//将该行的第 1 个数变成 1,即将等式两边同时除以这个数
for(int i=n; i>=c; i--){
a[r][i] /= a[r][c];
}
// 把第 c 列除了第 1 个数全部消成 0
for(int i=r+1; i<n; i++){
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; // 无解
}
}
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; // 有唯一解
}