电脑桌面
添加小米粒文库到电脑桌面
安装后可以在桌面快捷访问

北航数值分析大作业二VIP专享VIP免费

北航数值分析大作业二_第1页
北航数值分析大作业二_第2页
北航数值分析大作业二_第3页
数值分析大作业(二)学院名称宇航学院专业名称航空宇航推进理论与工程学生姓名段毓学号SY16153062016年11月5日1算法设计方案首先将矩阵A进行拟上三角化,把矩阵A进行QR分解,计算出RQ。要得出矩阵A的全部特征值,首先对A进行QR的双步位移得出特征值。最后,采用列主元的高斯消元法求解特征向量。1.1A的拟上三角化因为对矩阵进行QR分解并不改变矩阵的结构,因此在进行QR分解前对矩阵A进行拟上三角化可以大大减少计算机的计算量,提高程序的运行效率。具体算法如下所示,记,并记的第r列至第n列的元素为。对于执行若全为零,则令,转5;否则转2。计算令。计算继续。1.2A的QR分解具体算法如下所示,记,并记,令对于执行1.若全为零,则令,转5;否则转2。2.计算3.令。4.计算5.继续。1.3A的特征值为了加快收敛速度,采用带双步移位的QR分解求解A的全部特征值,具体算法如下,使用矩阵的拟上三角化算法把矩阵A化为你上三角阵;给定经度水平和最大迭代次数L。记,令k=1,m=n。如果,则得到A的一个特征值,置m:=m-1(降阶),转(5),否则转(4)。如果,则得到矩阵A的两个特征值,为二阶子阵的两个特征根,置m:=m-1(降阶),转(5);否则转(6)。如果m=1,则得到A的一个特征根,转(9),如果m=0,转(9)。否则转(3)。如果k=L,则计算失败。否则转(7)。记,计算,转(3)。A的全部特征值以计算完毕,停止计算。其中,的分解与的计算用下列算法实现:记。对于执行若全为零,则令,转5;否则转2。计算令。计算继续。此算法执行完后,就得到。1.4A的特征向量记λ为矩阵A的实特征值,x为对应的特征向量。则Ax=λx,即(A-λI)x=0的解即为矩阵A的特征向量。因此对于特征向量的求解可采用列主元素的Gauss消元法。其具体算法如下,消元过程对于执行选行号,使。交换与所含的数值。对于计算回代过程最终得到的向量的即为对应于实特征值的特征向量。2.程序运行结果2.1拟上三角化后2.2QR分解Q阵R阵RQ阵2.3QR的双步位移求特征值和特征向量特征值特征向量3源程序#include#include#definen11#defineE1e-12#defineL2500doublefuzhi(doublea[10][10]);//给矩阵A赋值doublenishangsanjiaohua(doublea[10][10]);doubleQR(doublea[10][10]);doublesgn(doublex);voidgauss();doubleA[10][10];doubletzR[n],I[n][n];intnR;voidQRshuangbu();voidmain(){fuzhi(A);nishangsanjiaohua(A);QR(A);QRshuangbu();gauss();}doublefuzhi(doublea[10][10])//给矩阵A赋值{inti,j;for(i=1;i<11;i++){for(j=1;j<11;j++){if(i==j){a[i-1][j-1]=1.52*cos(i+1.2*j);}else{a[i-1][j-1]=sin(0.5*i+0.2*j);}}}}doublesgn(doublex)//定义sgn函数{if(x>0){return1;}if(x<=0){return-1;}}doublenishangsanjiaohua(doublea[10][10]){intr,i,j;doublepingfanghe;doubledrpingfang;doubledr,cr,hr;doubleur[10],pr[10],qr[10],tr,wr[10];for(r=1;r<9;r++){pingfanghe=0;drpingfang=0;for(i=1;i<11;i++){ur[i-1]=0;pr[i-1]=0;qr[i-1]=0;wr[i-1]=0;}tr=0;//置零for(i=r+2;i<11;i++){pingfanghe=pingfanghe+a[i-1][r-1]*a[i-1][r-1];}if(pingfanghe==0){continue;//(1)}else{for(i=r+1;i<11;i++){drpingfang=drpingfang+a[i-1][r-1]*a[i-1][r-1];}dr=sqrt(drpingfang);//(2.1)drcr=-sgn(a[r][r-1])*dr;//(2.2)crhr=cr*cr-cr*a[r][r-1];//(2.3)hrfor(i=r+1;i<11;i++){ur[i-1]=a[i-1][r-1];}ur[r]=ur[r]-cr;//(3)urfor(j=1;j<11;j++){for(i=1;i<11;i++){pr[j-1]=pr[j-1]+a[i-1][j-1]*ur[i-1];}pr[j-1]=pr[j-1]/hr;}//(4.1)prfor(i=1;i<11;i++){for(j=1;j<11;j++){qr[i-1]=qr[i-1]+a[i-1][j-1]*ur[j-1];}qr[i-1]=qr[i-1]/hr;}//(4.2)qrfor(i=1;i<11;i++){tr=tr+pr[i-1]*ur[i-1];}tr=tr/hr;for(i=1;i<11;i++){wr[i-1]=qr[i-1]-tr*ur[i-1];}for(i=1;i<11;i++){for(j=1;j<11;j++){a[i-1][j-1]=a[i-1][j-1]-wr[i-1]*ur[j-1]-ur[i-1]*pr[j-1];}}}}for(i=1;i<11;i++){for(j=1;j<11;j++){printf("A%d,%d=%.12e\t",i,j,a[i-1][j-1]);}printf("\n");}}doubleQR(doublea[10][10]){doubleQ[10][10];doubleRQ[10][10];inti,j;intt;intr;doublepingfanghe;do...

1、当您付费下载文档后,您只拥有了使用权限,并不意味着购买了版权,文档只能用于自身使用,不得用于其他商业用途(如 [转卖]进行直接盈利或[编辑后售卖]进行间接盈利)。
2、本站所有内容均由合作方或网友上传,本站不对文档的完整性、权威性及其观点立场正确性做任何保证或承诺!文档内容仅供研究参考,付费前请自行鉴别。
3、如文档内容存在违规,或者侵犯商业秘密、侵犯著作权等,请点击“违规举报”。

碎片内容

文章天下+ 关注
实名认证
内容提供者

各种文档应有尽有

确认删除?
VIP
微信客服
  • 扫码咨询
会员Q群
  • 会员专属群点击这里加入QQ群
客服邮箱
回到顶部