高斯消元&&luogu3389

                            

高斯消元(Gauss)

高斯消元和我们做二元一次方程组差不多

流程:

1.把系数和右边的值就是用二维数组存下来->转化成矩阵

   我们的目标是把这个矩阵装换成 上三角的形式

对角线系数全部为1,1下面都为0,为了下面的回带

1 4 2 3
0 1 7 9 
0 0 1 2
0 0 0 1

2.利用 加减消元和等式两边除以一个数,一列一列的进行消元

   顺便判断一下是否有解,对角线上系数不为0

3.求出上三角之后,我们倒着回代一下就可以求取解了

 

当选取主元的时候,由于是double类型,当对角线的系数太小时,此时用它做除数会带来误差扩散,使结果严重失真。所以我们在消元的过程中,如果出现主元相差较大,要选取最大数作为主元,并交换行列,(当然,消元完毕的上边不能考虑在内)

 ---参考数学一本通


 

 代码

 1 #include <iostream>
 2 #include <cmath>
 3 #include <cstdio>
 4 using namespace std;
 5 
 6 const double eps=1e-6;
 7 int n;
 8 double a[110][110];
 9 double ans[110];
10 
11 int main()
12 {
13     scanf("%d",&n);
14     for(int i=1; i<=n; ++i)
15           for(int j=1; j<=n+1; ++j)
16             scanf("%lf", &a[i][j]);
17             
18     for(int i=1; i<=n; ++i) {
19         int pivot=i;
20         for(int j=i+1; j<=n; ++j)//选取较大主元 
21             if(fabs(a[j][i]) > fabs(a[pivot][i])) pivot=j;
22         if(abs(a[pivot][i]) < eps) { //判断有无解,无穷解也当做无解 
23             printf("No Solution");
24             return 0;
25         }
26         if(pivot!=i) swap(a[i],a[pivot]);//直接交换 
27         double tmp=a[i][i];
28         for(int j=i; j<=n+1; ++j) {
29             a[i][j]/=tmp;//系数化为1 
30         }
31         for(int j=i+1;j<=n;j++) {//下面的化为0 
32             tmp=a[j][i];
33             for(int k=i;k<=n+1;k++) {
34                 a[j][k]-=a[i][k]*tmp;
35             }
36         }
37     }
38     ans[n]=a[n][n+1];
39     for(int i=n-1; i>=1; i--) {
40         ans[i]=a[i][n+1];
41         for(int j=i+1; j<=n; ++j)
42            ans[i]-=a[i][j]*ans[j];
43     }//回带 
44     for(int i=1;i<=n;++i)
45         printf("%.2lf\n",ans[i]);
46 }
posted @ 2018-06-24 16:32  ComplexPug  阅读(184)  评论(0编辑  收藏  举报