高斯消元法

网友投稿 1057 2022-08-27

高斯消元法

高斯消元法

高斯消元法主要有消元和回代两个过程,和其他求解线性方程组的方法如克莱姆法则,约当消除法相比运算量小,速度快。更适合在计算机内使用。

假设有n个方程,var个变量,时间:

for i=1 --> n-1

for j=i+1-->n

for k=j+1-->n+1

所以时间是O(n^3)

把方程组写成增广矩阵,进行行初等变换,系数矩阵形成上三角矩阵,然后回代逐一求解每一个未知数。

结果又分为:无解;无穷多的解(找出自由变量);浮点数解;唯一整数解等情况

高斯消元的模板:

#include#include#include#include#includeusing namespace std;const int MAXN=50;int a[MAXN][MAXN];//增广矩阵int x[MAXN];//解集bool free_x[MAXN];//不确定的变元int equ,var,free_num;void Debug(){ int i, j; for (i = 0; i < equ; i++) { for (j = 0; j < var + 1; j++) { cout << a[i][j] << " "; } cout << endl; } cout << endl;}int gcd(int a,int b){ return b?gcd(b,a%b):a;}int lcm(int a,int b){ return a/gcd(a,b)*b;}// 高斯消元法解方程组(free_num=-2表示有浮点数解,但无整数解,//-1表示无解,0表示唯一解,大于0表示无穷解,并返回自由变元的个数)//有equ个方程,var个变元。增广矩阵行数为equ,分别为0到equ-1,列数为var+1,分别为0到var.int Gauss(int equ,int var){ int i,j,k; int max_r; int col;//当前处理的列 int ta,tb; int LCM; int temp; int free_x_num; int free_index; col=0; for(k = 0;k < equ && col < var;k++,col++) {// 找到该col列元素绝对值最大的那行与第k行交换.(为了处理无解的判断) max_r=k; for(i=k+1;iabs(a[max_r][col])) max_r=i; } if(max_r!=k) { for(j=k;j= 0; i--) { free_x_num = 0; // 用于判断该行中的不确定的变元的个数,如果超过1个,则无法求解,它为不确定的变元. for (j = 0; j < var; j++) { if (a[i][j] != 0 && free_x[j]) free_x_num++, free_index = j; } if (free_x_num > 1) continue; temp = a[i][var]; for (j = 0; j < var; j++) { if (a[i][j] != 0 && j != free_index) temp -= a[i][j] * x[j]; } x[free_index] = temp / a[i][free_index]; free_x[free_index] = 0; } return var - k; // 自由变元有var - k个. } for (i = var - 1; i >= 0; i--) { temp = a[i][var]; for (j = i + 1; j < var; j++) { if (a[i][j] != 0) temp -= a[i][j] * x[j]; } if (temp % a[i][i] != 0) return -2; // 说明有浮点数解,但无整数解. x[i] = temp / a[i][i]; } return 0;}void in(){ memset(a, 0, sizeof(a)); memset(x,0,sizeof(x)); memset(free_x,1,sizeof(free_x)); int i,j; for (i = 0; i < equ; i++) { for (j = 0; j < var + 1; j++)scanf("%d", &a[i][j]); }}void out(){ int i,j; if (free_num == -1) printf("no solutions!\n"); else if (free_num == -2) printf("here are double but not Integer solutions!\n"); else if (free_num > 0) { printf("here are many solutions and var number is %d\n", free_num); for (i = 0; i < var; i++) { if (free_x[i]) printf("x%d is not certain\n", i + 1); else printf("x%d: %d\n", i + 1, x[i]); } } else for (i = 0; i < var; i++)printf("x%d: %d\n", i + 1, x[i]); printf("\n");}int main(void){ freopen("cin.txt", "r", stdin); int i, j; while (scanf("%d %d", &equ, &var) != EOF) { in(); free_num = Gauss(equ,var); out(); } return 0;}

测试用例:

3 3

2 1 -1 8

-3 -1 2 -11

-2 1 2 -3

3 3

3 2 0 0

1 1 0 -1

0 0 0 0

3 3

3 2 0 0

1 1 1 -1

0 1 1 -2

3 3

3 2 0 0

1 1 0 -1

1 1 0 0

对应的答案是:

2

3

-1

2

-3

任意数

1

-3/2

-1/2

程序的结果: x1: 2 x2: 3 x3: -1 here are many solutions and var number is 1 x1: 2 x2: -3 x3 is not certain here are double but not Integer solutions! no solutions! 中间的过程(以:

3 3 2 1 -1 8 -3 -1 2 -11 -2 1 2 -3 为例): -3 -1 2 -11 2 1 -1 8 -2 1 2 -3 -3 -1 2 -11 0 1 1 2 0 5 2 13 -3 -1 2 -11 0 5 2 13 0 1 1 2 -3 -1 2 -11 0 5 2 13 0 0 3 -3 -3 -1 2 -11 0 5 2 13 0 0 3 -3 -3 -1 2 -11 0 5 2 13 0 0 3 -3 x1: 2 x2: 3 x3: -1

高斯消元过程仅仅涉及行变换,所以得到的X[]解集顺序是不变的。

压缩的写法:

#include #include #include #include using namespace std;const int N=105;int gcd(int a,int b){ return b==0?a:gcd(b,a%b);}// a[][]:增广矩阵 n: 方程数 m:变量数void Gauss(int a[N][N],int n,int m,int &r,int &c){ for(r=0,c=0;rabs(a[maxi][c])) maxi=i; } if(maxi!=r){ for(int j=r;j=0;i--){ int dex=0; int cnt=0; for(int j=0;j1) continue; int t=a[i][m]; for(int j=0;j=0;i--){ int t=a[i][c]; for(int j=i+1;j>n>>m){ for(int i=0;i

浮点数型的高斯消元法:

#include #include #include #include #include using namespace std;const int N=105;const double eps=1e-7;int myCmp(double a,double b){ if(fabs(a-b)b)-(a=0) a-=m; while(myCmp(a,0)<0) a+=m; return a;}double gcd(double a,double b){ return myCmp(b,0)==0?a:gcd(b,myMod(a,b));}// a[][]:增广矩阵 n: 方程数 m:变量数void Gauss(double a[N][N],int n,int m,int &r,int &c){ for(r=0,c=0;r0) maxi=i; } if(maxi!=r){ for(int j=r;j=0;i--){ int dex=0; int cnt=0; for(int j=0;j1) continue; double t=a[i][m]; for(int j=0;j=0;i--){ double t=a[i][c]; for(int j=i+1;j>n>>m){ for(int i=0;i

用几个式子验证,结果让人满意。

用mathematica计算:

In[1]:= Solve[{0.001*x1 + 2 x2 == 1, 3*x1 + 5*x2 == 2}, {x1, x2}]Out[1]= {{x1 -> -0.166806, x2 -> 0.500083}}In[2]:= Solve[{0.0000000000000001*x1 + 20000*x2 == 10000, 300*x1 + 500*x2 == 2000}, {x1, x2}]Out[2]= {{x1 -> 5.83333, x2 -> 0.5}}In[3]:= Solve[{0.0000000000000000001*x1 + 20000*x2 == 10000, 30000000000000*x1 + 500*x2 == 2000}, {x1, x2}]Out[3]= {{x1 -> 5.83333*10^-11, x2 -> 0.5}}高斯运算结果:-0.166805671392827 0.5000834028356965.833333333333333 0.5000000000000000.000000000058333 0.500000000000000Process returned 0 (0x0) execution time : 0.375 sPress any key to continue.

版权声明:本文内容由网络用户投稿,版权归原作者所有,本站不拥有其著作权,亦不承担相应法律责任。如果您发现本站中有涉嫌抄袭或描述失实的内容,请联系我们jiasou666@gmail.com 处理,核实后本网站将在24小时内删除侵权内容。

上一篇:ACdream 1084 寒假安排(数论(k进制末位的0的数目+n!中v因子个数))
下一篇:十佳最受欢迎的编程语言 你擅长几个?
相关文章

 发表评论

暂时没有评论,来抢沙发吧~