高斯列主元消去法C++代码
一、思想原理:
算法:
第1步 输入系数矩阵A,右端项b,置k:=1;
第2步 对k=1,…,n-1进行如下操作:
(1)选列主元,确定第一行第一个元素,保证它在第一列中绝对值是最大的,即为主元素;
(2)若第二行第二例的元素比第三行第二例的绝对值小,则交换的这两行;以此类推每一行的,直至保证每一行的对应位置的元素皆为绝对值最大的元素即可;
(3)消元,将矩阵运算化为上三角形即可
第3步 回代
将线性方程组做成增广矩阵,对增广矩阵进行行变换。
在第列中,第行及以下的元素选取绝对值最大的元素,将该元素最大的行与第行交换,然后采用高斯消元法将新得到的消去第行以下的元素,一次进行直到。从而得到上三角矩阵。再对得到的上三角矩阵进行回代操作,即可以得到方程组的解。
二、C++程序代码
#include<iostream> #include<cmath> using namespace std; #define N 10 //定义矩阵的维数,可自定义 static double A[N][N]; //系数矩阵 static double b[N]; //右端项 static double Aug[N][N + 1]; //增广矩阵 //高斯列主元消去法解线性方程组 void gaussin_LZY(double A[N][N],double b[N],double x[N],int n) { int i,i1, j, k; double Aug[N][N + 1], maxele, Temp, M; for (i = 0; i < n; i++) { //构造增广矩阵Aug for (j = 0; j < n; j++) Aug[i][j] = A[i][j]; Aug[i][n] = b[i]; } for (i1 = 0; i1 < n - 1; i1++) { //列主元 maxele = fabs(Aug[i1][i1]); //取主元,fabs是取绝对值函数 k = i1; for (i = i1; i < n; i++) if (maxele < fabs(Aug[i1][i]))//比较i1列 k = i1; //换行 for (j = 0; j < n + 1; j++) { Temp = Aug[i1][j]; Aug[i1][j] = Aug[k][j]; Aug[k][j] = Temp; } //消元,以第i1行为工具行处理下一行的元素 for (k = i1 + 1; k < n; k++) { M = -Aug[k][i1] / Aug[i1][i1]; for (j = i1; j < n + 1; j++) Aug[k][j] = Aug[k][j] + M * Aug[i1][j]; } } //回代求解 double s; x[n - 1] = Aug[n - 1][n] / Aug[n - 1][n - 1]; for (i = n - 2; i >= 0; i--) { s = 0; for (j = i + 1; j < n; j++) s = s + Aug[i][j] * x[j]; x[i] = (Aug[i][n] - s) / Aug[i][i]; }三、实例验证结果
//Aug=[1 1 1 1 10;-1 2 -3 1 -2;3 -3 6 -2 7;-4 5 2 -3 0] //x[i]={1 2 3 4}
四、MATLAB程序验证
function x=magauss2(A,b,flag) if nargin<3,flag=0;end n=length(b); for k=1:(n-1) %选主元 [ap,p]=max(abs(A(k:n,k))); p=p+k-1; if p>k t=A(k,:);A(k,:)=A(p,:);A(p,:)=t; t=b(k);b(k)=b(p);b(p)=t; end m=A(k+1:n,k)/A(k,k); %消元 A(k+1:n,k+1:n)=A(k+1:n,k+1:n)-m*A(k,k+1:n); b(k+1:n)=b(k+1:n)-m*b(k); A(k+1:n,k)=zeros(n-k,1); if flag~=0,Ab=[A,b],end end x=zeros(n,1); %回代 x(n)=b(n)/A(n,n); for k=n-1:-1:1 x(k)=(b(k)-A(k,k+1:n)*x(k+1:n))/A(k,k); end验证结果:结果一致,代码正确