普通矩阵特征值的QR算法.doc

普通矩阵特征值的QR算法 摘 要 求矩阵的特征值有多种不同的办法，本文主要介绍用QR法求矩阵的特征值，QR法是目前求中等大小矩阵全部特征值的最有效方法之一，使用于求实矩阵或复矩阵的特征值，它和雅可比法类似，也是一种变换迭代法。关键词：QR分解 迭代序列 特征值 Matlab一 、QR方法的理论： 对任意一个非奇异矩阵（可逆矩阵）A，可以把它分解成一个正交阵Q和一个上三角阵的乘积，称为对矩阵A的QR分解，即A=QR。如果规定R的对角元取正实数，这种分解是唯一的。若A是奇异的，则A有零特征值。任取一个不等于A的特征值的实数，则A-I是非奇异的。只要求出A-I的特征值和特征向量就容易求出矩阵A的特征值和特征向量，所以假设A是非奇异的，不是一般性。 设A=A1 ，对A1 作QR分解，得A1 = Q1R1，交换该乘积的次序，得A2 = R1Q1=，由于Q1正交矩阵，A1到A2的变换为正交相似变换，于是A1和A2就有相同的特征值。一般的令A1=A，对k=1，2，3，. 这样，可得到一个迭代序列Ak，这就是QR方法的基本过程。2、 QR方法的实际计算步骤Householder变换：如果 v 给出为单位向量而 I 是单位矩阵，则描述上述线性变换的是 豪斯霍尔德矩阵 （v * 表示向量 v 的共轭转置）H=I -2VV*三、化一般矩阵为Hessenberg阵 称形如的矩阵为上海森堡（Hessenberg）阵。如果此对角线元(i=2,3,.n)全不为零，则该矩阵为不可约的上Hessenberg矩阵。用Householder变换将一般矩阵A相似变换为Hessenberg矩阵。首先，选取Householder矩阵，使得经相似变换后的矩阵的第一列中有尽可能多的零元素。为此，应取为如下形式 其中为n-1阶Householder矩阵。于是有 只要取使得，就会使得变换后的矩阵的第一列出现n-2个零元。同理，可构造如下列形式Householder矩阵使得如此进行n-2次，可以构造n-2个Householder矩阵，使得。其中H为上Hessenberg矩阵。特别地，当A为实对称矩阵，则经过上述正交变换后，H变为三对角阵。用Householder方法对矩阵A作正交相似变换，使A相似于上Hessenberg阵，算法如下：（1） 计算（2） 计算（3） 计算四、上Hessenberg阵的QR分解对上Hessenberg阵只需要将其次对角线上的元素约化为零，用Given变换比用Householder变换更节省计算量。1、 平面旋转阵（Givens变换阵）定义 n阶方阵 称为平面旋转阵，或称为Givens变换阵。平面旋转阵的性质：（1） 平面旋转阵是非对称的正交阵。（2） 也是平面旋转阵。（3）平面旋转阵的作用：（1） 将向量的第j个分量约化为零。左乘向量只改变X的第i个分量和第j个分量。若令，有调整，可将约化为零。令，得所以，取，于是（2） 将向量的第个分量到第n个分量约化为零。（3） 用对矩阵A作变换得到的结论左乘A只改变A的第i，j行。右乘A只改变A的第i，j列。只改变A的第i，j行和第i，j列。2、用Givens变换对上Hessenberg阵作QR分解对上Hessenberg阵通常用n-1个Givens变换阵可将它化成上三角矩阵，从而得到B的QR分解式。具体步骤为： 设（否则进行下一步）取旋转矩阵则其中，设（否则进行下一步），再取旋转矩阵则其中。假设上述过程已进行了k-1步，有设，取其中，于是因此，最多做n-1次旋转变换，即得因为均为正交矩阵，故其中仍为正交矩阵。可算出完成这一过程的运算量约为4,比一般矩阵的QR分解的运算量少一个数量级。可证明仍是上Hessenberg阵，于是可按上述步骤一直迭代下去，这样的得到的QR方法的运算量比基本QR方法大为减少。具体实例现在我们就用上述所说的QR方法求解矩阵的全部特征值。第一步：将A化成上Hessenberg阵，取于是H即为与A相似的上Hessenberg阵。将H进行QR分解，记取于是再取于是 第一次迭代得重复上述过程，迭代11次得 从得出的结果我们可以看出其结果还是很接近精确值的，其实我们完全可用Matlab实现上述计算。用海森伯格QR基本算法求矩阵全部特征值Matlab得到运算结果为：源程序如下：function l = hessqrtz(A,M)%海森伯格QR基本算法求矩阵全部特征值%已知矩阵:A%迭代步数:M%求得的矩阵特征值:lA=5 -3 2;6 -4 4;4 -4 5;M=11;A = hess(A);for(i=1:M) q,r=qr(A); A = r*q; l = diag(A);End而用QR基本算法求矩阵全部特征值的运行结果如下：其Matlab源程序如下：function l = qrtz(A,M)%QR基本算法求矩阵全部特征值%已知矩阵:A%迭代步数:M%求得的矩阵特征值:lA=5 -3 2;6 -4 4;4 -4 5;M=11;for(i=1:M) %M为迭代步数 q,r=qr(A); A = r*q; l = diag(A);end 由以上两种运算结果的对比可以看