数值分析第二次大作业.docxVIP

题目分析使用带双步位移的QR分解法求矩阵的全部特使用带双步位移的QR分解法求矩阵的全部特征值，并对其中的每一个实特征值求相应的特征向量。已知： (i,j=1,2,……,10)算法设计按照题目给出的矩阵定义对矩阵A赋初值：对应的函数为initA()；为了减少求特征值和特征向量过程中的计算量，在对矩阵进行QR分解前先进行拟上三角化：对应的函数为nssj()；对拟上三角化后的矩阵A使用带双步位移的QR分解法逐次迭代（最大迭代次数L=500），逐个求出其特征值，对应的函数为tezhengzhi()；这个过程中包含着两个子程序：QR()和qmk()，分别用来对矩阵Mk进行QR分解并得到Ak+1和计算mk的值；使用带原点平移的反幂法求出其对应的特征向量，对应的函数为：fmifa()，这个过程中求解线性方程时用到列主元的高斯消元法，对应的函数为：gauss(double );根据数值分析课本的相关知识，步骤3中带双步位移的QR分解法的流程图如下：本次作业所使用的编译环境为：Visual c++6.0程序源代码#includestdio.h#includemath.hdouble A[10][10],rr[10][10],qq[10][10],rq[10][10],uk[10][10];double xy[10],x[10]={0};//反幂法中double tzr[10]={0},tzi[10]={0};//定义矩阵的特征值数组,r实部、i虚部//对矩阵A进行初始化，赋值voidinitA(){inti,j;for(i = 1; i = 10; i++) { for(j= 1; j= 10;j++) A[i-1][j-1]=sin(0.5*i+0.2*j); }for(i= 1; i= 10;i++)A[i- 1][i- 1]=cos(i+1.2*i)*1.5;}//定义函数对矩阵A进行拟上三角化voidnssj(){intr,i,j;intsgn(double a);double d,c,h,t,tr,u[10],p[10],q[10],w[10];double sum=0,tmp=0; //开始循环for (r=0;r8;r++) {t= 0;for(i= r+2;i10;i++) t= t+(A[i][r]==0);if (t==8- r) continue;else{sum=0;for (i=r+1;i10;i++)sum=sum+ A[i][r]*A[i][r]; d= sqrt(sum); c= -1*sgn(A[r+1][r])*d; h= c*c- c*A[r+1][r]; //step (3)for (i=0;i10;i++)u[i]= 0;for (i=r+2;i10;i++)u[i]= A[i][r];u[r+1]= A[r+1][r]- c; //step (4)for(i= 0;i 10; i++) {tmp= 0;sum= 0;for(j= 0;j 10; j++) {tmp= tmp+A[j][i]* u[j];sum= sum+A[i][j]* u[j]; }p[i]= tmp/h; q[i]= sum/h; }sum=0;for(i= 0;i 10; i++) sum= sum+ p[i]*u[i];tr= sum/h;for(i= 0;i 10; i++)w[i]= q[i]- tr*u[i];for(i= 0;i 10;i++) {for(j= 0;j 10;j++) { A[i][j]= A[i][j]- w[i]*u[j]- u[i]*p[j]; } }} } //以上的A[][]存在一定的误差，消除误差for(i=0; i10; i++) {for(j= 0; j 10; j++) if (-1.0e-12A[i][j] A[i][j] 1.0e-12) A[i][j] = 0; }}//sgn子程序intsgn(double a) {if(a = 0) return 1;else return -1;}//计算A的特征值的函数（课本P63 11步）inttezhengzhi(){int k=1,m=9,n=9;int