第3章:线性方程组求解 第3章:线性方程组求解代码汇编 问题:求Ax-b的解,A是M阶可逆方阵 约定:算法中用到的是M×N增广矩阵,N=M+1; 变量:i,j,k等为整型变量,x,y,z为实型变量; 1.把任意M阶线性方程组化为上三角形方程 组的c语言代码 void ToTriangleo I int i,j, k; float x. for(k=0; k<M; k++) 【x=A[k[k] forge=k:j<N: j++)A[k][]/=x; for(i=k+1; i<M; i ++ 【x=A[i][k for(=k;j<N;j++)A[][]-= x*A[k][]
第 3 章:线性方程组求解 1/5 第 3 章:线性方程组求解代码汇编 问题:求 Ax=b 的解,A 是 M 阶可逆方阵; 约定:算法中用到的是 M×N 增广矩阵,N=M+1; 变量:i,j,k 等为整型变量,x,y,z 为实型变量; 1. 把任意 M 阶线性方程组化为上三角形方程 组的 C 语言代码: void ToTriangle() { int i,j,k; float x; for(k=0;k<M;k++) { x = A[k][k]; for(j=k;j<N;j++) A[k][j] /= x; for(i=k+1;i<M;i++) { x=A[i][k]; for(j=k;j<N;j++)A[i][j] -= x*A[k][j]; }
第3章:线性方程组求解 2/5 return. 2.求解上三角形线性方程组c语言代码: for(k=M-1; k>0; k-) for(i=0; i<k;i ++) A[][N-1]-=A[i][k]*A[kIN-1 A[i]k]=0.0; 3约当消去法C语言代码: for(k=0 K<M; K++) [X=A[K[K: for(j=k<N;j++)A[K][]/=X for(i=0; i<M; i ++ 讦f(i==K) continue; X=AOKI
第 3 章:线性方程组求解 2/5 } return; } 2. 求解上三角形线性方程组 c 语言代码: for(k=M-1;k>0;k--) for(i=0;i<k;i++) { A[i][N-1] -= A[i][k]*A[k][N-1]; A[i][k]=0.0; } 3.约当消去法 C 语言代码: for(K=0;K<M;K++) { x=A[K][K]; for(j=K;j<N;j++) A[K][j] /= x; for(i=0;i<M;i++) { if(i==K) continue; x=A[i][K];
第3章:线性方程组求解 3/5 for(j=Kj<Nj++)A[[-=X*AK[j; 4利用约当消去法求矩阵的逆 对于MxM可逆矩阵A我们可以构造矩阵T=(AD,那么用约当消 去法把T的左边M列化为单位矩阵后其效果相当于用A1左乘 T从而得到A1(AD)=(A)所以T的后面的M列就是我们所要 求的A1 4实现选主元的源代码 j=K x=Math.abs(A[KIkd) for(=K+1;i<M++) i if(Math. abs(au[k)<x) continue x=Math.abs(Au][k]) 」=
第 3 章:线性方程组求解 3/5 for(j=K;j<N;j++)A[i][j] -=x*A[K][j]; } } 4.利用约当消去法求矩阵的逆 对于 M×M 可逆矩阵A,我们可以构造矩阵 T=(A I),那么用约当消 去法把 T 的左边 M 列化为单位矩阵后,其效果相当于用 A-1左乘 T,从而得到 A-1 (A I)=(I A-1 ),所以 T 的后面的 M 列就是我们所要 求的 A-1 4 实现选主元的源代码 j=K; x=Math.abs(A[K][K]); for(i=K+1;i<M;i++) { if(Math.abs(A[i][K])<x) continue; x=Math.abs(A[i][K]); j=i; }
第3章:线性方程组求解 //say(j="+j+",X="+A[K]) fori=k;i<N; i++) IX=A[KJO: AK的=A[j[ ADIO=x; 提示:把 Math. abs0改写为fab0即成为C语言代码。Say0是 用于调试的一个小脚本程序。 追过程c语言代码: U[o]=D[B[0] V[D]=C[o]/B[0] for(=1;<Ni++) X=B[-A[v[i-1] U[=D]-Ui-1]^A而/x
第 3 章:线性方程组求解 4/5 //say("j="+j+", x="+A[j][K]); for(i=K;i<N;i++) { x=A[K][i]; A[K][i]=A[j][i]; A[j][i]=x; } 提示:把 Math.abs()改写为 fabs()即成为 C 语言代码。say()是 用于调试的一个小脚本程序。 追过程 c 语言代码: U[0]=D[0]/B[0]; V[0]=C[0]/B[0]; for(i=1;i<N;i++) { x=B[i]-A[i]*V[i-1]; U[i]=(D[i]- U[i-1]*A[i])/x;
第3章:线性方程组求解 V0=C[/X //赶过程c语言代码: XN-1=UN-1} for(i=N-2;i>=0;-)x[=U[]-v[i]*x[i+1]
第 3 章:线性方程组求解 5/5 V[i]=C[i]/x; } //赶过程 c 语言代码: X[N-1]=U[N-1]; for(i=N-2;i>=0;i--)X[i]=U[i]-V[i]*X[i+1];