对于高斯消元法求解线性方程组,
我的理解就类似于我们在做数学题时的加减消元法,
只是把它写成一个通用的程序运算过程
对于一个线性方程组,我们从左往右每次将一列对应的行以下的元通过加减消元消去,
每个元的系数最终组成一个上三角矩阵,再倒序回带,求出答案
为了保证程序的可操作性,我们每次要将用来消去下面的元的数化为1,
再将下面的行每个元的系数同时减去主行的系数*扩大的倍数,
这时倍数即为该行要消去的元的系数
建议看一下《数学一本通》的内容,介绍的比较浅显
寻找主元:
double的除法操作是有一些误差的,我们在操作时,必须每次找一个开头最大的主元消去,这样可以减小精度误差
回代过程:
ans[i]=m[i][n+1]-m[i+1][i]*ans[i+1]-m[i+2][i]*ans[i+2]-m[i+3][i]*ans[i+3]……-m[n][i]*ans[n];(ans[i]/系数1=ans[i])
这实际上和我们加减消元是有很大相似之处的
模板题https://www.luogu.org/problemnew/show/P3389
代码:
1 #include<iostream> 2 #include<cstdio> 3 #include<cmath> 4 using namespace std; 5 const int MAXN = 1010; 6 int n; 7 double a[MAXN][MAXN],ans[MAXN]; 8 int main() 9 { 10 scanf("%d",&n); 11 for(int i=1;i<=n;i++) 12 for(int j=1;j<=n+1;j++) 13 scanf("%lf",&a[i][j]); 14 for(int i=1;i<=n;i++) 15 { 16 int now=i; 17 for(int j=i+1;j<=n;j++) 18 if(fabs(a[now][i])<fabs(a[j][i])) now=j; //找主元 19 if(now!=i) swap(a[i],a[now]); //换到当前行 20 double d=a[i][i]; 21 if(d==0){ //对角线上有0,则会出现有元无法被消去,无法得到唯一解 22 puts("No Solution"); 23 return 0; 24 } 25 for(int j=i;j<=n+1;j++) a[i][j]/=d; //主元化为1 26 for(int j=i+1;j<=n;j++) 27 for(int k=i+1;k<=n+1;k++) 28 a[j][k]-=a[j][i]*a[i][k]; //下方每行消去同列的元,该行也同时进行变换 29 for(int j=i+1;j<=n;j++) 30 a[j][i]=0; 31 } 32 for(int i=n;i>=1;i--) 33 { 34 ans[i]=a[i][n+1]; 35 for(int j=i+1;j<=n;j++) 36 ans[i]-=a[i][j]*ans[j]; 37 } 38 for(int i=1;i<=n;i++) 39 printf("%.2lf ",ans[i]+1e-9); 40 return 0; 41 }