隐式QR法求实矩阵地全部特征值matlab实现.docVIP

实用标准文案 精彩文档 隐式QR法求实矩阵的全部特征值matlab实现 要求：用matlab编写通用子程序，利用隐式QR法求实矩阵的全部特征值和特征向量。 思想：隐式QR法实质上就是将一个矩阵 Schur化，之后求解特征值就比较方便。而隐式QR法还需要用到household变换，以及上hessenberg变换。最后使用QR迭代，达到Schur化的结果。 步骤： 将矩阵A上hessenberg化（算法6.4.1），送而得到一个上hessenberg形矩阵H； 可约性判定，也就是判断次对角线元素是否非零，如果次对角线元素非零，则不可约。 Schur化，也就是通过QR迭代，将矩阵H变化成为某些次对角线元素变成0，同时还要满足，这些元素之间间隔最大为1，那么，所得到的最重的矩阵H就是一个Schur形矩阵。 假如两个等于0的次对角线元素间隔为0，那么该元素的上面一个元素，也就是H的对角线上的元素，即为其中一个特征值；假如两个等于0的次对角线元素间隔为1，那么在这两个元素之间就形成了一个2*2的矩阵，可以求解一个一元二次方程来得到两个共轭的特征值。 实验代码：详见附录2 实验结果：（代码相见附录2） （i）设矩阵A如下： 求x=0.9, 1.0, 1.1时的特征值和特征向量。 X=0.9： r是特征值，V是特征向量矩阵。 X=1： r是特征值，V是特征向量矩阵。 X=1.1： r是特征值，V是特征向量矩阵。 （ii） 求的所有根。 附录2 隐式QR迭代： 主程序： function [r,V]=SchurQR(A) %向量r用来储存特征值 %Hessenberg分解: [m,m]=size(A); for k=1:m-2 [v,b]=house(A(k+1:m,k)); H1= eye(m-k)-b*v*v; H2=eye(m); for i=k+1:m for j=k+1:m H2(i,j)=H1(i-k,j-k); end end if k==1; H=H2; else H=H*H2; end A(k+1:m,k:m)=H1*A(k+1:m,k:m); A(1:m,k+1:m)= A(1:m,k+1:m)*H1; end u=10e-5; for i=2:m; if abs( A(i,i-1))=(abs(A(i,i))+ abs(A(i-1,i-1)))*u; A(i,i-1)=0; end end %QR迭代： H22=A; x=Ifreducible(H22); while x==1 H22=Francis(H22); x=Ifreducible(H22); end [r,V]=EigValue(H22); 子程序1： function [r,V]=EigValue(A)%计算A的特征值，特征向量 [n,n]=size(A); r=zeros(1,n); y=zeros(1,n-1);%y用来储存次对角线元素 for i=1:n-1 y(i)=A(i+1,i); end m=0; for i=1:n-1 if abs(y(i)-0)1e-5 m=m+1; end end if m==0 x=1; else z=zeros(1,m);%z用来储存值为0的y向量的角标。 j=1; i=1; while(in) if abs(y(i)-0)1e-5 z(j)=i; j=j+1; end i=i+1; end end if z(1)==2%次对角线第一个等于0的元素的位置不同，需要2分类讨论 p=[1,A(1,1)+A(2,2),A(1,1)*A(2,2)-A(1,2)*A(2,1)] r(1:2)=roots(p);%求2*2矩阵的特征值 j=1; while jm if z(j+1)-z(j)==1 r(z(j+1))=A(z(j+1),z(j+1)); end if(z(j+1)-z(j)==2) p=[1,-(A(z(j+1)-1,z(j+1)-1)+A(z(j+1),z(j+1))),A(z(j+1)-1,z(j+1)-1)*A(z(j+1),z(j+1))-A(z(j+1)-1,z(j+1))*A(z(j+1),z(j+1)-1)]; r((z(j+