欢迎来到天天文库
浏览记录
ID:58854706
大小:323.50 KB
页数:18页
时间:2020-09-23
《数值分析报告第二大题.doc》由会员上传分享,免费在线阅读,更多相关内容在教育资源-天天文库。
1、学号:DY1305姓名:指导老师:《数值分析》计算实习报告第二大题一、题目要求已知10*10阶的矩阵A,使用带双步位移的QR分解法求其全部特征值,及每个实特征值对应的特征向量。此外还有拟上三角化产生的矩阵,及其进行QR分解结束后所得的矩阵Q、R以及乘积矩阵RQ。二、算法设计方案先用课本第61到62页的算法将矩阵进行拟上三角化,得到拟上三角化矩阵。再用课本第63到64页的带双步位移的QR算法求矩阵的全部特征值,同时可以得到矩阵A进行QR分解后的矩阵Q、R以及乘积矩阵RQ。由于矩阵A与矩阵相似,则他们的特征值相同,即可得矩阵A的全部特征值。对于各实特征值对应的特征向量,根据,利用实
2、特征值选主元Gauss消元法求解对应的特征向量。三、运算结果及分析以下计算结果截图。图1(a)对A进行拟上三角化得到的矩阵A(n-1)图1(b)对A进行拟上三角化得到的矩阵A(n-1)图2对A进行QR得到的矩阵Q图3对A进行QR得到的矩阵R图4对A进行QR得到的乘积矩阵RQ图5矩阵A的全部特征值图6各实特征值对应的特征向量四、源程序#include#include#include#include#include#defineEpsilon1e-12//精度#defineN11#d
3、efinen10#defineL1000//迭代次数doublea[N][N];doublelambda[N][2];doublefuzhi();//赋值函数voidNSSJH();//拟上三角化函数voidQR();//QR分解函数voidDoubleQR();//带双步位移的QR分解函数voidTZXL();//求取实特征根对应特征向量的函数voidmain(){voidNSSJH();voidQR();DoubleQR();TZXL();}doublefuzhi()//赋值函数,存储行列均为1-10{inti,j;for(i=1;i4、j++){if(i==j){a[i][j]=1.5*cos(i+1.2*j);}elsea[i][j]=sin(0.5*i+0.2*j);}}returna[i][j];}voidNSSJH()//拟上三角化函数{inti,j,r;doublecr,dr,hr,tr,temp;doublepr[N]={0},qr[N]={0},ur[N],wr[N]={0};fuzhi();for(r=1;r<=n-2;r++){dr=0;cr=0;hr=0;tr=0;temp=0;for(i=r+1;i<=n;i++){temp=temp+a[i][r]*a[i][r];}dr=sqrt(t5、emp);//求解drif(a[r+1][r]>0)cr=-dr;elsecr=dr;//求解crhr=cr*cr-cr*a[r+1][r];//求解hrfor(i=1;i<=r;i++){ur[i]=0;}ur[r+1]=a[r+1][r]-cr;for(i=r+2;i<=n;i++)//生成ur[N]{ur[i]=a[i][r];}for(i=1;i<=n;i++){temp=0;for(pr[i]=0.0,j=1;j<=n;j++)//求解pr[N]{pr[i]=pr[i]+a[j][i]*ur[j]/hr;}}for(i=1;i<=n;i++)//求解qr[N]{tem6、p=0;for(qr[i]=0.0,j=1;j<=n;j++){qr[i]=qr[i]+a[i][j]*ur[j]/hr;}}for(i=1;i<=n;i++)//求解tr{tr=tr+ur[i]*pr[i]/hr;}for(i=1;i<=n;i++)//求解wr[N]{wr[i]=qr[i]-ur[i]*tr;}for(i=1;i<=n;i++)//求解A(r+1){for(j=1;j<=n;j++){a[i][j]=a[i][j]-wr[i]*ur[j]-ur[i]*pr[j];}}}cout<<"******对A进行拟上三角化得到的矩阵A(n-1)******"<7、l;for(i=1;i<=n;i++){for(j=1;j<=n;j++){cout<<"A["<
4、j++){if(i==j){a[i][j]=1.5*cos(i+1.2*j);}elsea[i][j]=sin(0.5*i+0.2*j);}}returna[i][j];}voidNSSJH()//拟上三角化函数{inti,j,r;doublecr,dr,hr,tr,temp;doublepr[N]={0},qr[N]={0},ur[N],wr[N]={0};fuzhi();for(r=1;r<=n-2;r++){dr=0;cr=0;hr=0;tr=0;temp=0;for(i=r+1;i<=n;i++){temp=temp+a[i][r]*a[i][r];}dr=sqrt(t
5、emp);//求解drif(a[r+1][r]>0)cr=-dr;elsecr=dr;//求解crhr=cr*cr-cr*a[r+1][r];//求解hrfor(i=1;i<=r;i++){ur[i]=0;}ur[r+1]=a[r+1][r]-cr;for(i=r+2;i<=n;i++)//生成ur[N]{ur[i]=a[i][r];}for(i=1;i<=n;i++){temp=0;for(pr[i]=0.0,j=1;j<=n;j++)//求解pr[N]{pr[i]=pr[i]+a[j][i]*ur[j]/hr;}}for(i=1;i<=n;i++)//求解qr[N]{tem
6、p=0;for(qr[i]=0.0,j=1;j<=n;j++){qr[i]=qr[i]+a[i][j]*ur[j]/hr;}}for(i=1;i<=n;i++)//求解tr{tr=tr+ur[i]*pr[i]/hr;}for(i=1;i<=n;i++)//求解wr[N]{wr[i]=qr[i]-ur[i]*tr;}for(i=1;i<=n;i++)//求解A(r+1){for(j=1;j<=n;j++){a[i][j]=a[i][j]-wr[i]*ur[j]-ur[i]*pr[j];}}}cout<<"******对A进行拟上三角化得到的矩阵A(n-1)******"<7、l;for(i=1;i<=n;i++){for(j=1;j<=n;j++){cout<<"A["<
7、l;for(i=1;i<=n;i++){for(j=1;j<=n;j++){cout<<"A["<
此文档下载收益归作者所有