北航数值分析第二次大作业--QR分解.docxVIP

《数值分析A》计算实习题目二姓 名学 号联系方式班 级指导教师2012年10月算法设计方案整个程序主要分为四个函数，主函数，拟上三角化函数，QR分解函数以及使用双步位移求解矩阵特征值、特征向量的函数。因为在最后一个函数中也存在QR分解，所以我没有采用参考书上把矩阵M进行的QR分解与矩阵Ak的迭代合并的方法，而是在该函数中调用了QR分解函数，这样增强了代码的复用性，减少了程序长度；但由于时间关系，对阵中方法的运算速度没有进行深入研究。1.为了减少QR分解法应用时的迭代次数，首先对给定矩阵进行拟上三角化处理。2.对经过拟上三角化处理的矩阵进行QR分解。3.注意到计算特征值与特征向量的过程首先要应用前面两个函数，于是在拟上三角化矩阵的基础上对QR分解函数进行了调用。计算过程中，没有采用goto语句，而是根据流程图采用其他循环方式完成了设计，通过对迭代过程的合并，简化了程序的循环次数，最后在计算特征向量的时候采用了列主元高斯消去法。二、源程序代码#includestdio.h#includemath.h#includestring.hinti,j,k,l,m;//定义外部变量double d,h,b,c,t,s;double A[10][10],AA[10][10],R[10][10],Q[10][10],RQ[10][10]; double X[10][10],Y[10][10],Qt[10][10],M[10][10];double U[10],P[10],T[10],W[10],Re[10]={0},Im[10]={0};double epsilon=1e-12;void main(){voidQuasiuppertriangular(double A[][10]);voidQRdecomposition(double A[][10]);voidDoublestepsQR(double A[][10]);inti,j;for(i=0;i10;i++){for(j=0;j10;j++) { A[i][j]=sin(0.5*(i+1)+0.2*(j+1)); Q[i][j]=0; AA[i][j]=A[i][j]; } A[i][i]=1.5*cos(2.2*(i+1)); AA[i][i]=A[i][i];}Quasiuppertriangular(A); //调用拟上三角化函数printf( \n A经过拟上三角化矩阵为：\n\n);for(i=0;i10;i++) //输出拟上三角化矩阵{for(j=0;j10;j++) { printf(%.12e ,A[i][j]);//输出拟上三角化矩阵 }printf(\n\n);}QRdecomposition(A); //调用QR分解函数printf( 进行QR分解后，R矩阵为：\n\n); //输出R矩阵for(i=0;i10;i++) {for(j=0;j10;j++){ printf(%.12e ,R[i][j]);}printf(\n\n);}printf( Q矩阵为:\n\n);//输出Q矩阵for(i=0;i10;i++){for(j=0;j10;j++){ printf(%.12e ,Q[i][j]);}printf(\n\n);}printf( RQ矩阵为：\n\n);//输出RQ矩阵for(i=0;i10;i++){for(j=0;j10;j++){ printf(%.12e ,RQ[i][j]);}printf(\n\n);}DoublestepsQR(A);//调用双步位移函数printf( \n\n 特征值实部依次为：\n\n);//输出特征值实部for(j=0;j10;j++) {printf(%.12e ,Re[j]);}printf(\n\n 特征值虚部依次为：\n\n ); //输出特征值虚部for(j=0;j10;j++) {printf(%.12e ,Im[j]);}//按行输出特征向量printf( \n\n 按行输出实特征根相应特征向量为：\n\n);for(i=0;i10;i++){if(i==1||i==2||i==5||i==6) {continue; }for(j=0;j10;j++){ printf(%.12e ,X[i][j]);}printf(\n\n);}getchar();}//拟上三角化函数voidQuasiuppertriangular(double