标签:
高斯消元可能这名字听着挺高大上……但其实……没错他就是挺高大上的……他可以用来解线性方程组(不知道线性方程组是什么的自己去百度吧)!!!至于什么逆矩阵之类的我还没有研究。。。这里我先介绍一下主元高斯消元法
x-2y+3z=6(1)
4x-5y+6z=12(2)
7x-8y+10z=21(3)
这是一个三元一次方程组,先来手动解一下这个方程会对后面的方法有很大帮助!!!
我们看到这个方程,一定很希望让他变成如下的形式
x+0y+0z=a
0x+y+0z=b
0x+0y+z=c
所以我就要变形啦!!!先来(2)-(1)*4,然后再来一发(3)-(1)*7,这样除了(1)中的x都被消掉了。同理,我们可以得到唯一含有y的式子和唯一含有z的式子
x-2y+3z=6
3y-6z=-12(4)
6y-11z=-21(5)
现在我们的target是只让(4)中含有未知数y,消去其他方程中的y,首先把(4)/3
然后就得到了y-2z=-4(6)
最后我们用(1)-(6)*(-2)再来(5)-(6)*6
x-z=-2(7)
y-2z=-4
z=3(8)
然后随便算一算就知道(x=1,y=2,z=3)的方程的解了。这样即使未知数增加,也可以用同样的方法来求
下面我给出代码(其实我在学的时候觉得这个代码是不太好理解的…………
1 typedef vector<double> vec; 2 typedef vector<vec> mat; 3 4 vec gauss_jordan(const mat &A,const vec &b){ 5 int n=A.size(); 6 mat B(n,vec(n+1)); 7 for(int i=0;i<n;i++){ 8 for(int j=0;j<n;j++)B[i][j]=A[i][j]; 9 } 10 for(int i=0;i<n;i++)B[i][n]=b[i]; 11 //把b和A合并一发,这样比较好处理 12 for(int i=0;i<n;i++){ 13 int pivot=i; 14 for(int j=i;j<n;j++){ 15 if(abs(B[j][i])>abs(B[pivot][i]))pivot=j; 16 } 17 swap(B[i],B[pivot]); 18 //选择要消去的未知数系数的绝对值尽量大的方程可以减小误差 19 for(int j=i+1;j<=n;j++)B[i][j]/=B[i][i]; 20 //把正在处理的式子未知数系数变为1 21 for(int j=0;j<n;j++){ 22 if(i!=j){ 23 //第j个式子中消去第i个未知数 24 for(int k=i+1;k<=n;k++)B[j][k]-=B[j][i]*B[i][k]; 25 } 26 } 27 } 28 vec x(n); 29 for(int i=0;i<n;i++){x[i]=B[i][n];} 30 //最后存放在数组最右边的就是答案 31 return x; 32 }
其实这个完全可以用二维数组来实现……
)
标签:
原文地址:http://www.cnblogs.com/543Studio/p/5159748.html