雅可比迭代法与矩阵的特征值.docx

1、雅可比迭代法与矩阵的特征值实验五矩阵的lu分解法，雅可比迭代法班级： 学号：实验五 矩阵的LU分解法，雅可比迭代一、目的与要求：熟悉求解线性方程组的有关理论和方法；会编制列主元消去法、LU 分解法、雅可比及高斯塞德尔迭代法德程序；通过实际计算，进一步了解各种方法的优缺点，选择合适的数值方法。二、实验容：会编制列主元消去法、LU 分解法、雅可比及高斯塞德尔迭代法德程序，进一步了解各种方法的优缺点。三、程序与实例列主元高斯消去法算法：将方程用增广矩阵Ab=(表示1)消元过程对k=1,2,n-1选主元，找使得=如果，则矩阵A奇异，程序结束；否则执行。如果，则交换第k行与第行对应元素位置， j=k,n

2、+1消元，对i=k+1, ,n计算对j=l+1, ,n+1计算2)回代过程若，则矩阵A奇异，程序结束；否则执行。；对i=n-1, ,2,1,计算程序与实例程序设计如下：#include #include using namespace std;void disp(double* p,int row,int col) for(int i=0;irow;i+) for(int j=0;jcol;j+) coutpij ; coutendl; void disp(double* q,int n) cout=endl; for(int i=0;in;i+) coutXi+1=qiendl; cout=

3、endl;void input(double* p,int row,int col) for(int i=0;irow;i+) cout输入第i+1行:; for(int j=0;jpij; int findMax(double* p,int start,int end) int max=start; for(int i=start;iabs(pmaxstart) max=i; return max;void swapRow(double* p,int one,int other,int col) double temp=0; for(int i=0;icol;i+) temp=ponei;

4、ponei=potheri; potheri=temp; bool dispel(double* p,int row,int col) for(int i=0;irow;i+) int flag=findMax(p,i,row); /找列主元行号 if(pflagi=0) return false; swapRow(p,i,flag,col); /交换行 for(int j=i+1;jrow;j+) double elem=pji/pii; /消元因子 pji=0; for(int k=i+1;kcol;k+) pjk-=(elem*pik); return true;double sumRo

5、w(double* p,double* q,int row,int col) double sum=0; for(int i=0;i=0;i-) qi=(picol-1-sumRow(p,q,i,col)/pii; int main() coutn; double *p=new double* n; for(int i=0;in;i+) pi=new double n+1; input(p,n,n+1); if(!dispel(p,n,n+1) cout奇异endl; return 0;double* q=new doublen; for(int i=0;in;i+) qi=0; back(p

6、,n,n+1,q); disp(q,n); delete q; for(int i=0;in;i+) delete pi; delete p; 1. 用列主元消去法解方程例2 解方程组计算结果如下B=-1.461954C= 1.458125D=-6.004824E=-2.209018F= 14.719421矩阵直接三角分解法算法：将方程组Ax=b 中的A分解为A=LU，其中L为单位下三角矩阵，U为上三角矩阵，则方程组Ax=b化为解2个方程组Ly=b，Ux=y，具体算法如下：对j=1,2,3,n计算对i=2,3,n计算对k=1,2,3,n:a. 对j=k,k+1,n计算b. 对i=k+1,k+2

7、,n计算，对k=2,3,n计算,对k=n-1,n-2,2,1计算注：由于计算u的公式于计算y的公式形式上一样，故可直接对增广矩阵Ab=施行算法，此时U的第n+1列元素即为y。程序与实例例3 求解方程组Ax=bA=, b=结果为X0= 3.000001X1=-2.000001X2= 1.000000X3= 5.000000#include using namespace std;double* newMatrix(int row,int col) double *p=new double* row; /行 for(int i=0;irow;i+) /列 pi=new double col; re

8、turn p;void delMatrix(double* p,int row) for(int i=0;irow;i+) delete pi; delete p; void inputMatrix(double* p,int row,int col) for(int i=0;irow;i+) cout输入第i+1行:; for(int j=0;jpij; void dispMatrix(double* p,int row,int col) for(int i=0;irow;i+) for(int j=0;jcol;j+) coutpij ; coutendl; void dispVector

9、(double* q,int n) cout=endl; for(int i=0;in;i+) coutXi+1=qiendl; cout=endl;void initMatrix(double* p,int row,int col) for(int i=0;irow;i+) for(int j=0;j=col?col:row); double sum=0; for(int i=0;imin;i+) sum+=(Uicol*Lrowi); return sum;void resolve(double* A,double* L,double* U,int row,int col) for(int

10、 i=0;icol;i+) U0i=A0i; /把A的第一行给U的第一行 L00=A00; for(int i=1;irow;i+) /填充L的第一列 Li0=Ai0/A00; for(int i=1;irow;i+) for(int j=1;jcol;j+) if(i=j) Uij=Aij-sum(L,U,i,j); else Lij=(Aij-sum(L,U,i,j)/Ujj; for(int i=0;irow;i+) Lii=1;double sumRowY(double* L,double* y,int row) double sum=0; for(int i=0;irow;i+) s

11、um+=(Lrowi*yi); return sum;void backY(double* L,double* b,double* y,int row) for(int i=0;irow;i+) yi=0; /初始化y向量全为零 for(int i=0;irow;i+) yi=bi-sumRowY(L,y,i);double sumRowX(double* U,double* x,int row,int col) double sum=0; for(int i=row+1;icol;i+) sum+=(Urowi*xi); return sum;void backX(double* U,dou

12、ble* y,double* x,int row) for(int i=0;i=0;i-) xi=(yi-sumRowX(U,x,i,row)/Uii;int main() coutn; cout=endl; cout开始录入方阵数据.endl; double* A=newMatrix(n,n); /开辟矩阵A inputMatrix(A,n,n); /录入数据到A中 cout录入方阵数据完毕.endl; cout=endlendl; cout开始录入列阵.endl; cout输入列阵:; double* b=new doublen; /开辟向量b for(int i=0;ibi; /录入数据

13、到b中 cout录入列阵数据完毕.endl; double* L=newMatrix(n,n); /开辟方阵L initMatrix(L,n,n); /初始化L全为零 double* U=newMatrix(n,n); /开辟方阵U initMatrix(U,n,n); /初始化U resolve(A,L,U,n,n); /分解A为L与U double* y=new doublen; /开辟向量y backY(L,b,y,n); /回代求出y double* x=new doublen; /开辟向量X backX(U,y,x,n); /回代求出X dispVector(x,n); /释放空间

14、delMatrix(A,n); delMatrix(L,n); delMatrix(U,n); delete b; delete y; delete x;迭代法雅可比迭代法算法：设方程组Ax=b系数矩阵的对角线元素,M为迭代次数容许的最大值,为容许误差。取初始向量x=，令k=0。对i=1,2,n 计算如果，则输出，结束；否则执行。如果kM,则不收敛，终止程序；否则,转。程序与实例例4 用雅可比迭代法解方程组结果为迭代次数为20X0= 1.000000X1= 2.000000X2=-1.000000#include #include using namespace std;double* new

15、Matrix(int row,int col) double *p=new double* row; /行 for(int i=0;irow;i+) /列 pi=new double col; return p;void delMatrix(double* p,int row) for(int i=0;irow;i+) delete pi; delete p; void inputMatrix(double* p,int row,int col) for(int i=0;irow;i+) cout输入第i+1行:; for(int j=0;jpij; void dispMatrix(doubl

16、e* p,int row,int col) for(int i=0;irow;i+) for(int j=0;jcol;j+) coutpij ; coutendl; void dispVector(double* q,int n) for(int i=0;in;i+) coutXi+1=qit; coutendl;double sumRow(double* A,double* x1,int row,int col) double sum=0; for(int i=0;icol;i+) if(row=i) continue; sum+=(Arowi*x1i); return sum;void

17、iterat(double* A,double* b,double* x1,double* x2,int row) for(int i=0;irow;i+) x2i=(bi-sumRow(A,x1,i,row)/Aii;bool blow_error(double* x1,double* x2,int row,double e) double sum=0; for(int i=0;irow;i+) sum+=abs(x2i-x1i); if(sume) return true; else return false;int main() coutn; cout=endl; cout开始录入方阵数

18、据.endl; double* A=newMatrix(n,n); /开辟矩阵A inputMatrix(A,n,n); /录入数据到A中 cout录入方阵数据完毕.endl; cout=endlendl; cout开始录入列阵.endl; cout输入列阵:; double* b=new doublen; /开辟向量b for(int i=0;ibi; /录入数据到b中 cout录入列阵数据完毕.endl; cout=endl; coutM; coute; cout=endl; double* x1=new doublen; /开辟x1向量 double* x2=new doublen; /

19、开辟x2向量 cout输入初始向量:; for(int i=0;ix1i; cout=endl; int k=0; /迭代计数器 for(;) iterat(A,b,x1,x2,n); /迭代一次 if(blow_error(x1,x2,n,e) dispVector(x2,n); cout迭代次数为:k=M) cout不收敛endl; return 0; else k+; for(int i=0;in;i+) x1i=x2i; /释放空间 delMatrix(A,n); delete b; delete x1; delete x2;高斯-塞尔德迭代法算法:设方程组Ax=b的系数矩阵的对角线元

20、素,M为迭代次数容许的最大值,为容许误差取初始向量,令k=0。对i=1,2,n计算如果，则输出,结束；否则执行。如果则不收敛，终止程序；否则,转。程序与实例例5 用高斯-塞尔德迭代法解方程组结果为X0=3.000000X1=2.000000X2=1.000000#include #include #include #include #define N 100main() int i; float *x; float c12=8.0,-3.0,2.0,20.0, 4.0,11.0,-1.0,33.0, 6.0,3.0,12.0,36.0; float *GauseSeidel(float *,i

21、nt); x=GauseSeidel(c,3); for (i=0;i=2;i+) printf(x%d=%fn,i,xi); getch();float *GauseSeidel(float *a,int n) int i,j,nu=0; float *x,dx,d; x=(float *)malloc(n*sizeof(float); for (i=0;i=n-1;i+) xi=0.0; do for (i=0;i=n-1;i+) d=0.0; for (j=0;j=N) printf(fold numbern); nu+; while (fabs(dx)1e-6); return x;例

22、6 用雅可比迭代法解方程组迭代4次得解,若用高斯-塞尔德迭代法则发散。#includevoid main (void) int k,n; double x3=7,2,5; for(k=0;k5;k+) double a,b; a=x0; b=x1; x0=(7-2.0*x1+2*x2)/1; x1=(2-a-x2)/1; x2=(5-2*a-2*b)/1; for(n=0;n3;n+) printf(x%d=%8.6fn,n,xn); 用高斯-塞尔德迭代法解方程组迭代84次得解，若用雅克比迭代法则发散。#include#includevoid LOOP(float a1010,float b10,float x10,int);void main(void) float a1010