高斯消元法

Posted 2020-11-06 wxyww

高斯消元法

可以用于求解线性方程组，即n元1次方程组。利用矩阵，大致思路与普通解方程方法类似。只是更具一般性。将系数与右侧的常数存成一个矩阵，然后每次用第i行消去下面每行的第i个系数，最后就会得到一个一元方程，然后从后到前依次代回即可。

　　然后就是精度的问题，因为计算机中没有分数，所以只能用double来存，但是double的精度也是有限的，所以要想办法维护精度。

　　第一种优化办法就是在用第i行去消第i个系数时，从下面所有行中选择第i个系数的绝对值最大的一个与当前行交换，这样可以使每次消元时差距尽量大，然后就可以保持精度了。

　　第二种优化就是在进行消元时，不是用一个变量来储存两行首元素的商，而是每次直接计算。这样就必须防止第一行提前被修改，所以倒着进行修改。

代码：

 1 #include<cstdio>
 2 #include<cstring>
 3 #include<iostream>
 4 #include<cmath>
 5 using namespace std;
 6 const int N=110,M=1010;
 7 int n,m;
 8 double a[N][N],tmp[N][N],ans[N];
 9 inline void work()
10 {
11     for(int i=1;i<=n;++i)
12     {
13         int kk=i;
14         for(int j=i+1;j<=n;++j)
15             if(abs(a[kk][i])<abs(a[j][i]))
16                 kk=j;
17         if(kk!=i)
18         for(int j=1;j<=n+1;++j)
19             swap(a[kk][j],a[i][j]);
20         for(int j=i+1;j<=n;++j)
21         {
22             for(int k=n+1;k>=i;--k)
23                 a[j][k]=a[j][k]-a[j][i]/a[i][i]*a[i][k];
24         }
25     }
26     
27     for(int i=n;i>=1;--i)
28     {
29         for(int j=n;j>i;--j)
30         {
31             a[i][j]*=ans[j];
32             a[i][n+1]-=a[i][j];
33         }
34         ans[i]=a[i][n+1]/a[i][i];
35     }
36 }
37 int main()
38 {
39     scanf("%d%d",&n,&m);
40     for(int i=1;i<=n;++i)
41         for(int j=1;j<=n;++j)
42             scanf("%lf",&tmp[i][j]);
43     for(int k=1;k<=m;++k)
44     {
45         memset(ans,0,sizeof(ans));
46         memset(a,0,sizeof(a));
47         for(int i=1;i<=n;++i)
48             for(int j=1;j<=n;++j)
49                 a[i][j]=tmp[i][j];
50         for(int j=1;j<=n;++j)
51             scanf("%lf",&a[j][n+1]);
52         work();
53         for(int i=1;i<=n;++i)
54             printf("%.9llf ",ans[i]);
55         printf("\n");
56     }
57     return 0;
58 }