数值分析（B）大作业（二）1、算法设计：① 矩阵的拟上三角化：对实矩阵A进行相似变换化为拟上三角矩阵A（nT）,其变换矩阵采用Householde矩阵，变换过程如下：若a（r）=0（i=r,2,,n），则Hirr否则，s=（0,,0,a（r）,a（r））T，rr+1,rn,rc=-sgn（a（r））llsII，,a(r))，nrrr+1,r.r..-c,a(r),r+1,rrr,2,ru=s-ce=（0,,0,a（r）H=I-2u%t炖*，rrrr2A（r,1）=HA（r）Hrr=HHA(i)HH，n-211n-2当r=n—2时，得A（n-1）=・a（n€2）h=-n-2令P=HHH又H是对称正交矩阵，于是An-1=PtAP成立，因而An-1与A相12n-2r似•② 矩阵的QR分解：矩阵的QR分解过程与拟上三角化过程相似，在这里不再重复其原理③ 求全部特征值矩阵拟上三角化后利用带双步位移的QR方法，采用书本Page63页具体算法实现为了使编程方便，并体会goto语句使用的灵活性，程序中的跳转均使用gotoLoop*实现④ 求A的相应于实特征值的特征向量求实特征值对应的特征向量，即是求解线性方程组SPA）广°,（=I，n）。
因此，为得到J全部实特正值对应的特征向量，解线性方程组的过程要循环n次（n为矩阵阶数）线性方程组的求解采用列主元素Gauss消去法•••#include#include#defineERR1.0e-12//误差限#defineN10//矩阵行列数#defineL1.0e5//最大迭代次数doubleA[N][N]={0};voidInit_A()//初始化矩阵{inti,j;for(i=0;i=0)return1;elsereturn-1;}voidOn_To_The_Triangle()//矩阵拟上三角化{inti,j,r,flag=0;doublecr,dr,hr,tr,temp;doubleur[N],pr[N],qr[N],wr[N];for(r=1;r<=N-2;r++){flag=0;for(i=r+2;i<=N;i++)if(A[i-1][r-1]!=0){flag=1;break;}if(0==flag)continue;dr=0;for(i=r+1;i<=N;i++)dr+=A[i-1][r-1]*A[i-1][r-1];dr=sqrt(dr);if(0==A[r][r-1])cr=dr;elsecr=-Sgn(A[r][r-1])*dr;hr=cr*cr-cr*A[r][r-1];for(i=1;i<=r;i++)ur[i-1]=0;ur[r]=A[r][r-1]-cr;for(i=r+2;i<=N;i++)ur[i-1]=A[i-1][r-1];for(i=1;i<=N;i++){temp=0;for(j=1;j<=N;j++)temp+=A[j-1][i-1]*ur[j-1];pr[i-1]=temp/hr;}for(i=1;i<=N;i++){temp=0;for(j=1;j<=N;j++)temp+=A[i-1][j-1]*ur[j-1];qr[i-1]=temp/hr;}temp=0;for(i=1;i<=N;i++){temp+=pr[i-1]*ur[i-1];tr=temp/hr;}for(i=1;i<=N;i++){wr[i-1]=qr[i-1]-tr*ur[i-1];}for(i=1;i<=N;i++)for(j=1;j<=N;j++)A[i-1][j-1]=A[i-1][j-1]-wr[i-1]*ur[j-1]-ur[i-1]*pr[j-1];}}voidGet_Roots(doubleeigenvalue[][2],intm,doubless,doublett)//求一元二次方程的根{doublediscriminant=ss*ss-4*tt;//if(discriminant<0){*(*(eigenvalue+m-2))=0.5*ss;*(*(eigenvalue+m-2)+1)=0.5*sqrt(-discriminant);*(*(eigenvalue+m-1))=0.5*ss;*(*(eigenvalue+m-1)+1)=-0.5*sqrt(-discriminant);}else{*(*(eigenvalue+m-2))=0.5*(ss+sqrt(discriminant));*(*(eigenvalue+m-2)+1)=0;*(*(eigenvalue+m-1))=0.5*(ss-sqrt(discriminant));*(*(eigenvalue+m-1)+1)=0;}}voidGet_Mk(doublemk[][N],intm,doubless,doublett)〃获取Mk,用于带双步位移的QR分解{inti,j,k;for(i=0;i0)printf("+%8.4f",*(*(value+i)+1));elseif(*(*(value+i)+1)<0)printf("%8.4f",*(*(value+i)+1));printf("\n");}printf("\n");}intQR_With_Double_Step_Displacement(doubleeigenvalue[][2])//带双步位移QR分解求特征值{inti,j,k=1,m=N;doubles,t;doubleMk[N][N];for(i=0;ifabs(A[max][k]))max=i;if(k==max)return;else{for(i=k;i=0;k--){*(Eigenvector+k)=Temp_M[k][N];for(i=k;icnX-3.3B30JJhenb=1.577SiEigenuector=0.6217-0<1115-2.4838-1.3E569-3.81568.1173-1.2392-Q.68SEI2.6919WhenX=-1.48^0Eigenu巳ctor=Pl?.307824408190.4134-«.5?2OU.ld929-0.S374-B.34B3-0.3?8SUhenX=6.9356WhenX=-Eigpenu巳ctor=4.745S3.157917.299E-1.9803-31.87527.7946-10.042616.7B7613.1052UPen入=6.@565Eigenuector--S.1050-4i88S3-0.S788-^.6043-3.W4B715.7487-?.395H-7.1tJ97fZisjcnucctor"-0.43t7-0；9BG4Eicjcnvector—2.79241.5982-1.6679-12.25717.2412-5.398228.4101-12.1652-1.0724-1.0829-1.2&930.2G251.68292.1127PS:Forfurtherdetails,pleasecontactme:。