作用
高斯消元可以用来解线性方程组,也有一些别的用途。
算法流程
我们来看这样的一个方程组:
$ x+y+z=3$
$2x+3y+z= 6$
$4x+y+2z=7$
然后我们把每一项的系数提取出来,就得到如下的矩形 :
$1 \ 1 \ 1 \ 3 $
$ 2 \ 3 \ 1 \ 6 $
$ 4 \ 1 \ 2 \ 7 $
然后,我们就可以用矩阵的第一行乘$(a[i][1]\div a[1][1])$加到矩阵的第2行、第3行,使矩阵的2~3行的第1项变为0。如图:
$ 1 \ \ \ \ 1 \ \ \ \ 1 \ \ \ \ 3 $
$ 0 \ \ \ \ 1 \ -1 \ \ \ \ 0 $
$ 0 \ -3 \ -2 \ -5 $
然后我们再用第2行乘3并加到第三行,使第3行的第2项边为0。 如图:
$ 1 \ \ \ \ 1\ \ \ \ 1 \ \ \ \ 3 $
$ 0 \ \ \ \ 1 \ -1 \ \ \ \ 0 $
$ 0 \ \ \ \ 0 \ -5 \ -5 $
然后我们发现矩阵第3行只有一个未知数$z$,于是我们就可以直接求出$z$ 的值,然后第2行就只有一个未知数$y$了。求出$y$后第一行就只有一个未知数$x$。这样,我们就求出了这个方程组的解。
所以说,我们需要通过某些操作使得矩形变成只有一个倒三角上的数字可以不为0,然后依次求解每一个未知数的值。
注: 对于一个有n个未知数的方程组,若在矩阵处理完毕后前n列有至少一列上的数字全部为0,则该方程有多解,即这一行系数所对应的边量可以取任意值。
代码
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58
| #include<iostream> #include<cstring> #include<cstdio> #include<cmath> #include<algorithm> #define ll long long #define INF 0x7fffffff #define re register #define qwq printf("qwq\n");
using namespace std;
int read() { register int x = 0,f = 1;register char ch; ch = getchar(); while(ch > '9' || ch < '0'){if(ch == '-') f = -f;ch = getchar();} while(ch <= '9' && ch >= '0'){x = x * 10 + ch - 48;ch = getchar();} return x * f; }
double a[105][105],Max,x[105]; int n,flag;
void gauss(int n) { int k; for(int i = 1; i < n; i++) { k = i; Max = 0; for(int j = i + 1; j <= n; j++) if(Max < fabs(a[j][i])) Max = a[j][i], k = j; if(k != i) for(int j = i; j <= n + 1; j++) swap(a[k][j],a[i][j]); for(int j = i + 1; j <= n; j++) { double b = (a[j][i] / a[i][i]); for(int k = i; k <= n + 1; k++) a[j][k] = a[j][k] - a[i][k] * b; } } for(int i = n; i >= 1; i--) { if(a[i][i] == 0) {printf("No Solution\n"); flag = 1; return ;} for(int j = i + 1; j <= n; j++) a[i][n + 1] = a[i][n + 1] - a[i][j] * x[j]; x[i] = a[i][n + 1] / a[i][i]; } }
int main() { n = read(); for(int i = 1; i <= n; i++) for(int j = 1; j <= n + 1; j++) scanf("%lf",&a[i][j]); gauss(n); if(flag) return 0; for(int i = 1; i <= n; i++) printf("%.2f\n",x[i]); return 0; }
|