最近在看惯导的东西,然后想要用C++解惯导控制方程,然后就重头把C++解方程组这方面的知识回顾了一下,首先就是高斯列主元消去法,这个方法还算实用,这里以3*3的矩阵为例,里面注释很详细,各位小白可以参考参考hhh
//顺序高斯消去法求线性方程组的解 AX=B,A为系数矩阵,B为常向量矩阵,X为矩阵的解 //为了减小误差,采用列主元消去法,从第k列的a[k][k]及其以下的各元素中选取绝对值最大的元素,然后行变换将这两行交换位置 //对于非奇异矩阵并且消元过程中第K步的主对角线上的元素不为0,就可以用高斯消元法来求解 #include<iostream> const int n=3; void Gauss_elimination(double a[n][n],double b[n]) { double x[n];//解的存储数组 int i,j,k; double aa[n];//左上的主对角元素为0时系数矩阵的行相互交换 double bb;//左上的主对角元素为0时系数矩阵的行相互交换 double mm;//找出最大的a[m][k]行的a[m][k] int m;//找出最大的a[m][i]行的行m值 double c[n];//存储初等行变换的系数,用于行的相减 for(k=0;k<n-1;k++) //消元的整个过程如下,总共n-1次消元过程。 { mm=a[k][k];m=k;bb=0; for(int i=1;i<n-k;i++)//找出最大的绝对值a[m][k]行的行m值 { if(mm<abs(a[k+i][k])) {m=k+i;mm=abs(a[k+i][k]);} } if(mm==0) //此种情况表明a[k][k]以及下面全部是0,则方程有多解或者无解,然后结束函数 { std::cout<<"此方程无解或多解"<<std::endl; return; } bb=b[m]; b[m]=b[k]; b[k]=bb; for(int i=k;i<n;i++)//将最大的a[m][k]行移到最上面 { aa[i]=a[m][i]; a[m][i]=a[k][i]; a[k][i]=aa[i]; } for(i=k+1;i<n;i++)//求出第K次初等行变换的系数 c[i]=a[i][k]/a[k][k]; for(i=k+1;i<n;i++) //第K次的消元计算 { for(j=0;j<n;j++) { a[i][j]=a[i][j]-c[i]*a[k][j]; } b[i]=b[i]-c[i]*b[k]; } } std::cout<<"输出变换后的系数矩阵:"<<std::endl;//输出变换后的系数矩阵及常向量矩阵 for(int i=0;i<n;i++) { for(int j=0;j<n;j++) { std::cout<<a[i][j]<<" "; } std::cout<<"b["<<i+1<<"]="<<b[i]<<std::endl; } //先计算出最后一个未知数; x[n-1]=b[n-1]/a[n-1][n-1]; //求出每个未知数的值 for(i=n-2;i>=0;i--) { double sum=0; for(j=i+1;j<n;j++) { sum+=a[i][j]*x[j]; } x[i]=(b[i]-sum)/a[i][i]; } std::cout<<"此方程组的解为:"<<std::endl; for(int i=0;i<n;i++) { std::cout<<"x["<<i+1<<"]="<<x[i]<<std::endl; } }主函数如下:
int main() { double a[3][3]={1,1,1,0,4,-1,2,-2,1}; double b[3]={6,5,1}; Gauss_elimination(a,b); system("pause"); return 0; }运行结果如下: