数值分析QR方法求矩阵特征值与特征向量.docVIP

四.实验代码： function [H,B]=Hessenberg(A) n=length(A); B=eye(n); for k=1:n-2 X=zeros(n-k,1); H=eye(n); for i=1:n-k X(i)=A(i+k,k); end a=max(abs(X)); if a==0.0 break end X=X/a;c=X(1); b1=sqrt(sum(X.^2)); if X(1)=0 b1=-b1; end X(1)=X(1)-b1; b=b1^2-b1*c; H0=eye(n-k)-X*X/b; for i=1:n-k for j=1:n-k H(i+k,j+k)=H0(i,j); end end A=H*A*H; B=B*H;end H=A; 一. 实验题目： QR 方法求矩阵的特征和特征向量 二.设计目的： 学会利用镜面变换进行矩阵的QR分解及利用将幂法求 特征值和特征向量，熟悉Matlab编程环境。 三.设计原理： 利用镜像变换将A相似变换为Hessenberg B矩阵。记录变换矩阵。 运用Householder矩阵进行QR分解，QR方法为： B1=B B1=Q1R1 B2=R1Q1 .. . . Bm=QmRm Bm+1=RmQm Bm+1与Bm相似，从而特征值相等。 再利用原点位移的反幂法求B（或A）的特征向量。 反幂法用来计算矩阵按模最小的特征值及其特征向量，也可用来计算对应与一个给定近似特征值的特征向量。设A∈Rn×n 为非奇异矩阵，A的特征值依次记为 |λ1 |≥|λ2 |≥|λ3 |≥…≥|λn |,相应的特征向量为x1 ,x2 ,…,xn , 则A-1 的特征值为|1/λn |≥|1/λn-1 |≥…≥|1/λ1 | ,相应的特征向量为xn ,x n-1 ,…,x1 . 所以计算A的按模最小的特征值λn的问题就是计算A-1 的按模最大的特征值问题。对于A-1 应用幂法迭代(称为反幂法)，可求得矩阵A-1 的主特征值1/λn ，从而求得A的按模最小的特征值λn 。 反幂法迭代公式任取初始向量ν0 =μ0 ≠0，构造向量序列 (2.10)迭代向量νk 可以通过解方程组 Aνk =μk-1 求得.在反幂法中也可以用原点平移法来加速迭代过程或求其他特征值及特征向量。 如果矩阵(A-pI)-1 存在, 对其应用幂法，得反幂法的迭代公式 反幂法迭代公式中的νk 是通过解方程组 (A - pI)νk = μk-1 求得的。为节省工作量，可先将A - pI进行三角分解  (A - pI) = LU ,其中P为某个排列阵，于是求νk 相当于解两个三角形方程组 (2.12)实验表明，按下述方法选择μ0 较好：选μ0 使 Uν1 = L-1 μ0 = (1,1,…,1)T 用回代求解(2.12)即得ν1 ，然后再按公式(2.11)进行迭代。迭代公式：1. 分解计算 (A - pI) = LU , 且保存L,U信息。 2. 反幂法迭代 (1) 解 Uν1 =(1,1,…,1)T 求ν1  μ1 =ν1 / (ν1 )r ,  (2) k=2,3,… a) 解 Lyk =μk-1 ，求得yk 解 Uνk = yk ，求得νk b) μk =νk /(νk)r 步骤与实现：（1）利用镜像变换将A相似变换为Hessenberg B矩阵。记录变换矩阵 （2）利用QR方法求矩阵B的特征值。 （3）利用原点位移的反幂法求B（或A）的特征向量