QR迭代法求矩阵特征值
基于QR分解迭代求解方阵特征值和特征向量

基于QR分解迭代求解方阵特征值和特征向量一、特征值与特征向量求解的难点线性代数的知识告诉我们如果要求一个方阵的特征值,只需要求解如下的特征方程的根即可:但是在具体程序中,如何去求解一个高次的多项式方程的根本身就是一个难点,它的实现甚至要比求得特征值还要复杂。
因此,线性代数中这种用来手算的方法很明显不适合程序实现。
那么,我们该如何去编写程序求解特征值呢?二、QR迭代法求解特征值的原理经过搜集相关资料,已知在计算机编程领域求解特征值一般有两种比较常见的算法,一种是雅可比(Jacobi)迭代法,另一种就是我们要讲的QR迭代法,这两种算法本质上都是去不断迭代,通过逼近的方法去求得特征值。
雅可比方法有一个弊端,就是他只能求解实对称矩阵的特征值,适用面比较窄,但是QR分解迭代法可以求解任意矩阵的特征值,而且实现上更加接近我们线性代数课程的知识,因此下面介绍QR分解迭代法求解特征值的原理。
理论依据:任意一个非奇异矩阵(满秩的方阵)都可以分解为一个正交矩阵和一个上三角矩阵的乘积,且当对角元符号确定时,分解是唯一的。
QR分解是一种迭代方法,迭代格式如下:当基本上收敛到上三角矩阵的时候,迭代完成,此时主对角线元素就是特征值。
特别的:当是实对称矩阵时,是对角阵,就是其标准正交的特征向量矩阵,且有。
那么我们该如何去理解上述原理呢?下面进行简要的推导。
不妨令,对进行QR分解得到:令,很明显,对再进行QR分解:继续令,可以看到,同理,我们不断去进行QR分解,最终就会得到:因此,与拥有相同的特征值。
特别的,当时,有,将记作,那么和就分别是正交的特征向量矩阵和对角线元素为特征值的对角阵。
有了上述理论介绍,我们已经知道了如何去进行QR迭代求特征值,下面马上又有一个问题摆在了我们面前——如何用程序去做QR分解呢?三、QR分解的初等变换法QR分解有三种常用的算法:Gram-Schmidt正交化(线性代数讲授的方法)、初等变换法、Givens 变换法。
数值分析QR方法求矩阵特征值和特征向量

四.实验代码: function [H,B]=Hessenberg(A) n=length(A);B=eye(n);for k=1:n-2X=zeros(n-k,1);H=eye(n);for i=1:n-kX(i)=A(i+k,k);enda=max(abs(X));if a==0.0breakendX=X/a;c=X(1);b1=sqrt(sum(X.^2));if X(1)>=0b1=-b1;endX(1)=X(1)-b1;b=b1^2-b1*c;H0=eye(n-k)-X*X'/b;for i=1:n-kfor j=1:n-kH(i+k,j+k)=H0(i,j);end endA=H*A*H;B=B*H;endH=A;一.实验题目:QR方法求矩阵的特征和特征向量二.设计目的:学会利用镜面变换进行矩阵的QR分解及利用将幂法求特征值和特征向量,熟悉Matlab编程环境。
三.设计原理:利用镜像变换将A相似变换为Hessenberg B矩阵。
记录变换矩阵。
运用Householder矩阵进行QR分解,QR方法为:B1=BB1=Q1R1B2=R1Q1....Bm=QmRmBm+1=RmQmBm+1与Bm相似,从而特征值相等。
再利用原点位移的反幂法求B(或A)的特征向量。
反幂法用来计算矩阵按模最小的特征值及其特征向量,也可用来计算对应与一个给定近似特征值的特征向量。
设A∈R n×n为非奇异矩阵,A 的特征值依次记为|λ1|≥|λ2|≥|λ3|≥…≥|λn |,相应的特征向量为x1 ,x2,…,x n,则A-1的特征值为|1/λn|≥|1/λn-1|≥…≥|1/λ1 | ,相应的特征向量为x n ,x.所以计算A的按模最小的特征值λn的问题就是计算n-1,…,x1A-1的按模最大的特征值问题。
对于A-1应用幂法迭代(称为反幂法),可求得矩阵A-1的主特征值1/λn,从而求得A的按模最小的特征值λn。
用QR算法求矩阵的特征值

一.实验名称:用QR 算法求矩阵的特征值 二.实验目的:1、通过实验进一步熟悉掌握求矩阵特征值的QR 方法及原理。
2、 理解QR 方法的计算流程。
3、 能够编程实现QR 方法°A 和H 矩阵的全部特征值。
四. 实验要求:(1)根据QR 算法原理编写程序求矩阵A 及矩阵H 的全部特征值(要求误差<10* )。
(2)直接用MATLAB 的内部函数ei 吕求矩阵A 及矩阵H 的全部特征值,并与(1)的 结果比较。
五、QR 方法计算矩阵特征值的程序:functiontzl=diag(A);[namda, time, data_na]=qr_tz(A, tol) wuc ha=norm(t z 0-1 z1);if nargin==l;% 2 r三、实验内容:给左矩阵A =2 3 1,<1 i L(2 434 45 5 66、7H = 0 3 6 7 8 ,采用QR 方法计算0 0 2 8 9<0 0 0 1 0丿A=A1;tol=le-5;endwcha=l;time=O;while (wucha>tol)&(time<500)[q, r] =qr (A);Al=r*q;time=time+l;data_na(time, :)=tzl; endnamda^tzl; disp('特征值为')namdadisp('第一个特征在值')tzO=diag(Al);time迭代次数为time =nl=length(data_na);n2=(l :nD ,;templ=[n2, data_na];subplot (2, 2, 1:2)plot (date_na(:, 1))title('迭代次数为')gridsubplot (2, 2, 3)六、实验结果:» A=[6, 2, 1;2,3, 1;1,1, 1]; [namda,特征值为plot(data-na(:,2)) title('第二个特征值’) gridsubplot (2, 2, 4)plot(data-na(:,3)) title —第三个特征值') gridtime, data_na]=qr_tz (A, le-5);Q Figure 1File Edit View Insert Tools Desktop Window Help乂DdH^I第三个特征值Array» A=[6, 2, 1 ;2, 3, 1; 1,1, 1]; [V, D] =eig(A,' nobalance'),A二[2, 3, 4, 5, 6;4, 4, 5, 6, 7;0, 3, 6, 7, 8;0, 0, 2, & 9;0, 0, 0, 1, 0j; [namda, time, data_na]=qr_t z(A, le-5);特征值为namda =迭代次数为time =22»A= [2, 3, 4, 5,6;4, 4, 5, 6,7;0, 3, 6,7,8;0, 0, 2, 8,9;0, 0, 0,1, 0] ; [V, D]=eig(A,' nobalance1), V =0 00 00 0 0 0 0 0 0 00 0表1用两种方法求得矩阵A的全部特征值2H从图1和图2中可以看出在迭代前几次可能会有一些波动,但逐渐趋于平稳,并且收敛速度快,算法稳左。
特征值的快速求解技巧论文

特征值的快速求解技巧论文特征值问题是线性代数中一个重要的问题,涉及到矩阵的特征值和特征向量的求解。
求解特征值问题对于很多应用领域来说都是必不可少的,比如在机器学习、信号处理和图像处理中。
由于特征值问题的复杂性,研究人员一直在寻找快速求解特征值的方法。
本文将介绍一些特征值问题的快速求解技巧。
一、QR方法QR方法是一种常用的求解特征值问题的数值算法,它的基本思想是通过QR分解将一个矩阵转化为上三角矩阵,从而求解特征值。
QR方法可以分为隐式QR方法和显式QR方法。
隐式QR方法是通过不断迭代的方式求解特征值,在每次迭代中都会进行QR分解。
显式QR方法是直接计算出QR分解的结果,然后通过计算得到特征值。
二、幂方法幂方法是一种迭代法,主要用于求解矩阵的最大特征值及其对应的特征向量。
它的基本思想是选择一个初始向量,通过不断迭代矩阵的幂次,最终得到矩阵的最大特征值和对应的特征向量。
幂方法是一种简单而有效的求解特征值的方法,但它只能求解最大特征值及其对应的特征向量。
三、反迭代方法反迭代方法是对幂方法的改进,主要用于求解矩阵的特征值和特征向量。
它的基本思想是通过对幂方法的迭代结果应用逆迭代进行修正,从而得到更精确的特征值和特征向量。
反迭代方法通常比幂方法收敛更快,并且可以求解非对称矩阵的特征值和特征向量。
四、雅可比方法雅可比方法是一种常用的求解对称矩阵特征值的方法,它的基本思想是通过不断迭代地将矩阵中的非对角元素置为0,从而得到对称三对角矩阵。
然后,通过QR分解求解对称三对角矩阵的特征值,最终得到原矩阵的特征值。
雅可比方法可以保证得到的特征值是精确的,但是它的计算复杂度比较高,不适用于大规模矩阵的特征值求解。
五、奇异值分解方法奇异值分解方法是一种数值稳定的特征值求解方法,它可以求解任意矩阵的特征值和特征向量。
奇异值分解方法的基本思想是将矩阵分解为三个矩阵的乘积,其中一个矩阵是对角矩阵,其对角线上的元素即为矩阵的奇异值,另外两个矩阵分别是特征向量矩阵和其转置。
正交变换与QR迭代矩阵特征值计算

正交变换在QR迭代中的优势与局限性
优势
正交变换可以保证QR迭代的稳定性和收敛性,同时能够避免直接处理特征值问题,简化计算过程。
局限性
正交变换需要构造一个正交矩阵,计算量较大,对于大规模矩阵计算效率较低。此外,对于某些特殊 类型的矩阵,QR迭代可能不是最优的方法。
04
矩阵特征值计算
特征值的基本概念
稳定性好
相较于其他特征值计算方法,正交变换与QR迭代具有较 好的数值稳定性,能够有效地抵抗计算过程中的误差累积 ,从而得到更可靠的结果。
并行化潜力
正交变换与QR迭代中的各个步骤可以独立进行,这使得 这些方法具有并行化的潜力,能够利用高性能计算资源加 速特征值计算。
正交变换与QR迭代在矩阵特征值计算中的未来发展方向
收敛性条件包括迭代向量的范数逐渐减小以及特征值与迭代向量的夹角逐渐增大等。
03
正交变换在QR迭代中的应用
正交化过程在QR迭代中的重要性
保证迭代过程的稳定性
通过正交化过程,可以消除QR迭代中的数值误差,提高计算精度 和稳定性。
避免特征值问题
正交变换可以将原矩阵转化为上三角矩阵,从而避免直接处理特征 值问题,简化计算过程。
特征值
矩阵的特征值是矩阵的一个重要属性,它表示矩 阵对某些向量具有特定的缩放能力。
特征向量
与特征值对应的非零向量。
特征多项式
用于求解特征值的方程。
特征值的计算方法
代数法
01
通过求解特征多项式来找到特征值。
迭代法
02
通过迭代过程来逼近特征析矩阵的谱来找到特征值。
特征值计算中的问题与挑战
正交变换与QR迭代矩阵 特征值计算
• 正交变换概述 • QR迭代算法介绍 • 正交变换在QR迭代中的应用 • 矩阵特征值计算 • 正交变换与QR迭代在矩阵特征值计算
用qr方法求矩阵 的全部特征值例题

用qr方法求矩阵的全部特征值例题矩阵的特征值问题是矩阵理论中的重要问题之一,QR方法是一种常用的求解矩阵特征值的方法。
本文将通过一个具体的例题,介绍如何使用QR方法求矩阵的全部特征值。
一、问题描述给定一个$n\timesn$矩阵$A$,我们需要求出其全部特征值。
矩阵的特征值通常可以通过求解矩阵的特征多项式来得到。
对于实对称矩阵,我们可以通过对角化矩阵的方法来求解特征值。
但对于一般矩阵,我们需要使用其他方法,如QR方法。
二、QR方法原理QR方法是基于矩阵的QR分解原理,将原矩阵分解为一个正交矩阵和一个上三角矩阵的乘积。
通过这个分解,我们可以将原矩阵的特征多项式转化为一个简单的多项式,从而方便地求解特征值。
三、例题及解答【例题】给定一个$3\times 3$矩阵:$A=\begin{bmatrix}1&2&3\\0&-2&4\\0&-1&2\end{bmatrix}$要求求出该矩阵的全部特征值。
【解法】1. 将矩阵A进行QR分解,得到正交矩阵$Q$和上三角矩阵$R$:$A=QR$2. 计算$Q^TAQ$的特征多项式,并求出全部特征值。
3. 将上三角矩阵$R$代入特征多项式中,得到原矩阵A的特征值。
【代码实现】(使用MATLAB)```matlab% 定义矩阵AA = [1 2 3; 0 -2 4; 0 -1 2];% 进行QR分解[Q, R] = qr(A);% 计算Q^TAQ的特征多项式,并求出全部特征值[eigvals,~,~] = eig(Q^TAQ);eigvals = real(eigvals); % 取实部作为特征值% 将上三角矩阵R代入特征多项式中,得到原矩阵A的特征值eigenvalues = diag(R) ./ (diag(R)+eigvals);```【结果】经过以上步骤,我们可以得到原矩阵A的全部特征值为:$\lambda_1=2.75+0. 866i,\lambda_2=2.75-0.866i,\lambda_3=4$。
QR迭代法求矩阵特征值

本文用QR迭代法求解矩阵A的特征值:第一步先用豪斯荷尔德变换将矩阵A化为上海 森伯格矩阵AH ,第二步再对AH 进行QR迭代(使用吉文斯变换),当迭代满足精度要求时 输出n个特征值。
1
问题描述
对于给定的矩阵A(存在n个不相等的实特征值),采用QR方法求其n个特征值。本文中的矩 阵A由文件”gr9 009 00c rg.mm”给出。并以.mm格式描述。
2 AHi,i 2 + AHi +1,i
si = 于是,第i行和第i+1行的元素变换为:
AHi+1,i
2 + AH 2 AHi,i i+1,i
AHi,j = ci ∗ AHi,j + si ∗ AHi+1,j AHi+1,j = ci ∗ AHi+1,j − si ∗ AHi,j j = 1, 2, · · · · · · , n . 对i从1至n-1,重复以上过程,则等效实现QR分解,并得到上三角矩阵R。 在第二个过程中,R矩阵右乘上Q矩阵,即: AHnew = RQ = RG(1, 2; θ1 )T · · · · · · G(n − 2, n − 1; θn−1 )T G(n − 1, n; θn−1 )T 等效于依次对R矩阵的第i列和第i+1列进行列变换。变换后的元素为: Rj,i = ci ∗ Rj,i + si ∗ Rj,i+1 Rj,i+1 = ci ∗ Rj,i+1 − si ∗ Rj,i j = 1, 2, · · · · · · , n . 对i从1至n-1,重复以上过程,则等效地实现AHnew = RQ过程。 至此,完成了一次QR迭代。 本函数(AHnew=QR-iteration-once(AH,n))的流程图如图三所示。
QR方法QR方法是求任意矩阵的全部特征值的一种有效方法,它是

QR方法QR方法是求任意矩阵的全部特征值的一种有效方法,它是QR方法QR方法是求任意矩阵的全部特征值的一种有效方法,它是JACOBI方法的推广。
基本思想利用矩阵的QR分解,通过逆序相乘产生对原矩阵的一系列正交相似变换,使其变化为一个近似的上三角矩阵来求全部特征值。
这里QR分解是指将矩阵化为一个正交矩阵Q和一个上三角矩阵左乘的形式。
构造原理实对称矩阵可用正交相似变换将其化为对角形矩阵,但对非对称矩阵,一般用正交相似变换化不成对角矩阵,但SCHUR分解定理给我们一个有关这方面的结果。
定理3。
(实SCHUR分解定理)设矩阵A∈R n*n,则存在一个正交矩阵Q∈R n*n,使Q T AQ=其中每个B ii是1*1或2*2的小矩阵,若B ii为1*1的,其元素就是A的实特征值,否则B ii的特征值是A一对共轭复特征值。
此定理的证明可参阅文献[3]。
定理3指出了求矩阵A的全部特征值也可用正交相似变换的方法来做,正交相似变换的结果虽然不是对角矩阵,而是分块三角形矩阵,但它同样能很方便地求出全部特征值,有关一般矩阵的正交相似变换,我们不加证明地给出一个结论。
定理4。
设非奇异矩阵A∈R n*n,且有n个不同的特征值,记A=A(1)。
如果对整数k,有矩阵A(k)的QR分解为A(k)=Q k R k,则令A(k+1)=Q T k A(k)Q k,当k→∞时有A(k)本质上收敛于分块上三角形矩阵,这里“本质上收敛”指A(k)的主对角线上的元素或子块有确定的极限,其它元素或子块不管是否有极限。
此定理给出了求解一般矩阵全部特征值的方法。
由定理3,A(k+1)=(Q1Q2....Q k)T A(Q1Q2....Q k),令,则Q k也是正交矩阵,A(k+1)=Q T k A(k)Q k说明A(k+1)也是原矩阵A的正交相似变换,从而A(k+1)与A有相同的特征值,n任意,此外,由A(k)=Q k R k,则有Q T k A(k)= Q T k A(k)R k=R k,故有A(k+1)=Q k R k,这说明A(k+1)可直接交换Q k与R k的乘积顺序得到,于是可的如下QR算法。
qr迭代法求特征值特征向量

qr迭代法求特征值特征向量QR迭代法是求解特征值和特征向量的一种常用方法。
在这篇文章中,我们将详细介绍QR迭代法的原理、步骤和应用。
希望读者通过本文的阅读,能够对QR迭代法有更深入的理解。
QR迭代法是一种基于矩阵分解的数值方法,它可以用来求解一个实方阵或复方阵的特征值和特征向量。
QR迭代法的基本思想是通过不断地对矩阵进行正交相似变换,将矩阵转化为上三角矩阵或者对角矩阵,从而得到矩阵的特征值和特征向量。
QR迭代法的优点是收敛速度快,精度高,适用于大规模矩阵求解。
QR迭代法的步骤如下:1.将待求解的矩阵A分解为QR,其中Q是正交矩阵,R是上三角矩阵。
2.计算R的逆矩阵。
3.计算B=R的逆矩阵乘以A。
4.将B再次分解为QR。
5.重复步骤2-4,直到矩阵收敛为止。
在实际应用中,QR迭代法可以用来求解各种问题,例如线性方程组的求解、特征值分解、奇异值分解等。
下面我们分别介绍一下这些应用。
1.线性方程组的求解对于一个线性方程组Ax=b,可以通过QR分解和反向代入来求解。
具体步骤如下:1.将系数矩阵A进行QR分解,得到Q和R两个矩阵。
2.将b进行正交变换,得到新的向量y=Q^Tb。
3.利用R和y求解新的线性方程组Rx=y。
4.通过反向代入得到x的解。
2.特征值分解对于一个实对称矩阵A,可以通过QR迭代法来求解其特征值和特征向量。
具体步骤如下:1.将A进行QR分解,得到Q和R两个矩阵。
2.计算B=R的逆矩阵乘以Q^T,得到新的矩阵C=B^T乘以A 乘以B。
3.重复步骤1-2,直到C收敛为止。
4.将C的对角线元素作为A的特征值,B的列向量作为A的特征向量。
3.奇异值分解对于一个m*n的实矩阵A,可以通过QR迭代法来进行奇异值分解。
具体步骤如下:1.将A的转置矩阵AT与A进行乘积运算得到ATA。
2.对ATA进行QR分解,得到正交矩阵Q和上三角矩阵R。
3.计算V=Q。
4.对A进行乘积运算得到AT*A。
5.对AT*A进行QR分解,得到正交矩阵Q和上三角矩阵R。
矩阵特征值的求法举例

矩阵特征值的求法举例特征值是线性代数中一个重要的概念,它能够描述一个矩阵对应的线性变换的特性。
在实际应用中,我们经常需要计算一个矩阵的特征值。
本文将通过举例来讲解矩阵特征值的求法。
我们来介绍一下什么是特征值。
给定一个n×n的矩阵A,如果存在一个非零向量v,使得Av=λv,其中λ为常数,那么我们称λ为矩阵A的特征值,而v称为矩阵A对应于特征值λ的特征向量。
计算矩阵特征值的方法有很多,包括特征值分解、幂法、反幂法、QR方法等。
下面我们来逐一介绍这些方法,并通过具体的例子进行说明。
1. 特征值分解法特征值分解是指将一个矩阵分解成特征值和特征向量的乘积的形式,即A=QΛQ^-1,其中Q是特征向量组成的矩阵,Λ是对角矩阵,其对角线上的元素是矩阵A的特征值。
举例:假设有一个2×2的矩阵A=[4, 2; 1, 3],我们来计算其特征值。
首先我们要求解方程det(A-λI)=0,其中I是单位矩阵,λ是待求的特征值。
展开方程可得(4-λ)(3-λ)-2·1=0,解这个二次方程可得λ1=5,λ2=2。
2. 幂法幂法是一种迭代法,用于求解特征值模最大的特征值和对应的特征向量。
举例:假设有一个3×3的矩阵A=[1, 2, 3; 1, 3, 2; 3, 2, 1],我们来计算其特征值和特征向量。
首先我们随机选取一个初始向量x^(0),计算向量序列x^(k+1)=Ax^(k),迭代到收敛后,我们取得到的向量x^(k+1)的模最大的分量作为矩阵A的特征值模最大的特征向量。
然后,我们将这个特征向量归一化,即除以特征值模最大的分量,得到单位特征向量。
我们将单位特征向量与矩阵A相乘,可得到特征值l。
通过幂法计算可得矩阵A的特征值l≈2.863,以及对应的特征向量v≈[0.618, 0.618, 0.486]。
3. QR方法QR方法是一种迭代法,用于求解特征值。
举例:假设有一个5×5的矩阵A=[3, -1, 0, 0, 0; -1, 3, -1, 0, 0; 0, -1, 3, -1, 0; 0, 0, -1, 3, -1; 0, 0, 0, -1, 3],我们来计算其特征值。
qr迭代法求特征值特征向量

qr迭代法求特征值特征向量1. 介绍qr迭代法是一种数值计算方法,用于求解矩阵的特征值和特征向量。
特征值和特征向量是矩阵在线性代数中的重要概念,它们可以帮助我们理解矩阵的性质和行为。
qr迭代法是一种迭代算法,通过不断迭代矩阵的qr分解来逼近矩阵的特征值和特征向量。
2. qr分解qr分解是将一个矩阵分解为一个正交矩阵和一个上三角矩阵的乘积。
qr分解可以用于求解线性方程组、最小二乘问题等。
对于一个n阶矩阵A,它的qr分解可以表示为A=QR,其中Q是一个n阶正交矩阵,R是一个上三角矩阵。
3. qr迭代法的基本思想qr迭代法的基本思想是通过不断迭代矩阵的qr分解,使得矩阵收敛到一个上三角矩阵。
具体步骤如下:1.初始化一个n阶矩阵A和一个单位矩阵Q。
2.对矩阵A进行qr分解,得到A=QR。
3.计算新的矩阵A’ = RQ。
4.重复步骤2和步骤3,直到矩阵A’收敛到一个上三角矩阵。
5.上三角矩阵的对角线元素就是矩阵A的特征值的近似值。
4. qr迭代法的收敛性qr迭代法的收敛性是保证qr迭代法能够得到矩阵的特征值的一个重要性质。
qr迭代法的收敛性可以通过矩阵的谱半径来判断。
谱半径是矩阵所有特征值的绝对值的最大值。
如果矩阵的谱半径小于1,那么qr迭代法是收敛的。
5. qr迭代法的算法实现qr迭代法的算法可以通过以下步骤实现:1.初始化一个n阶矩阵A和一个单位矩阵Q。
2.设置一个收敛条件,比如矩阵的谱半径小于某个阈值。
3.while 矩阵的谱半径大于阈值:–对矩阵A进行qr分解,得到A=QR。
–计算新的矩阵A’ = RQ。
–更新矩阵A为A’。
4.上三角矩阵的对角线元素就是矩阵A的特征值的近似值。
6. qr迭代法的应用qr迭代法在科学计算和工程应用中有广泛的应用。
它可以用于求解大型稀疏矩阵的特征值和特征向量,求解线性方程组,计算矩阵的条件数等。
qr迭代法的收敛速度较快,精度较高,因此在实际应用中得到了广泛的应用。
7. 总结qr迭代法是一种求解矩阵特征值和特征向量的有效方法。
8.4 求矩阵全部特征值的QR方法

平面旋转阵Ri,j的性质: (1)RiT, j Ri , j I , Ri,1j RiT, j , 平面旋转阵是非对称的 正交阵。 T (2)Ri , j也是平面旋转阵。 (3)det(Ri , j )=1
(1)将向量x = x1 , x2 , ..., xn 的第j个分量约化为零。 T Ri , j 左乘向量x = x1 , x2 , ..., xn 只改变x的第i个分量和 第j个分量。 y1 x1 若令y Ri , j x,有 y R1,2 x 2 2 yi xi cos x j sin cos sin x1 x y j xi sin x j cos sin cos 2 y x k 1, ..., n; k i , j k , k x1 cos x2 sin 调整,可将y j约化为零。 x1 sin x2 cos xj 令y j 0,得 tan xj xi
用Householder阵 作正交相似变换
*
* *
二、化一般矩阵为上Hessenberg阵
称形如 h1n 1 h11 h12 h h2n 1 21 h22 h32 h33 H hnn 1 的矩阵为上海森堡(Hessenberg) 线元hii 1 ( i 2, 3, 约的上Hessenberg矩阵。 讨论用 Householder 变换将一般矩阵A相似变 换成Hessenberg阵 h1n h2 n h3 n hnn 阵。如果此对角
T a11 a2 H 1 于是有 H1 AH1 H 1a1 H1 A22 H 1 其中a1 (a21 , a31 , , an1 )T , a2 (a12 , a13 , , a1n )T ,
正交变换与QR迭代矩阵特征值计算

为Householder变换,或初等反射矩阵。
性质
(1)
对称:H
T w
Hw
(2)
正交:H
1 w
H
T w
Hw
(3)
对合:H
2 w
I
(4) 保模: Hwx 2
x 2
(5) det( Hw ) 1
3
Householder 变换
定理:设 x, y Rn, x y 且 ||x||2 = ||y||2,则存在 n 阶 Householder 变换 H,使得
R2
m
Rmm
拟上三角矩阵
其中 Rii 是一阶或二阶方阵。
若 Rii 是一阶方阵,则它就是 A 的特征值; 若 Rii 是二阶方阵,则其特征值为 A 的两个共轭复特征值。
13
QR 迭代
QR 迭代算法
计算矩阵的所有特征值和特征向量 计算过程
(1) 令 A1=A (2) 对 k = 1, 2, ... ,
y = Hx
证:取 w ( x y)
|| x y ||2
4
Householder 变换
定理:对任意的非零向量 x Rn,存在 Householder 变
换 H,使得
Hx = e1
其中 = sgn(x1)||x||2,
e1=
(1,
0,
...,
0)T
,H
I
1
uuT
u x e1
(1) 构造 H1 使得 H1 a1 = 1e1 ,令
1
a(2) 12
QR分解迭代求特征值——原生python实现(不使用numpy)

QR分解迭代求特征值——原⽣python实现(不使⽤numpy)QR分解:有很多⽅法可以进⾏QR迭代,本⽂使⽤的是Schmidt正交化⽅法具体证明请参考链接迭代格式实际在进⾏QR分解之前⼀般将矩阵化为上hessnberg矩阵(奈何这个过程⽐较难以理解,本⼈智商不够,就不做这⼀步了哈哈哈)迭代终⽌条件看了很多⽂章都是设置⼀个迭代次数,感觉有些不是很合理,本来想采⽤A(k+1)-A(k)的对⾓线元素的⼆范数来作为误差的,但是我有没有⼀些严格的证明,所以本⽂也采⽤⽐较⼤众化的思路,设置迭代次数。
Python实现1 M = [[2, 4, 2], [-1, 0, -4], [2, 2, 1]]23import copy4import math567class QR(object):89def__init__(self, data):10 self.M = data11 self.degree = len(data)1213def get_row(self, index):14 res = []15for i in range(self.degree):16 res.append(self.M[i][index])17return res1819def get_col(self, index):20 res = []21for i in range(self.degree):22 res.append(self.M[i][index])23return res2425 @staticmethod26def dot(m1, m2):27 res = 028for i in range(len(m1)):29 res += m1[i] * m2[i]30return res3132 @staticmethod33def list_multi(k, lt):34 res = []35for i in range(len(lt)):36 res.append(k * lt[i])37return res3839 @staticmethod40def one_item(x, yArr):41 res = [0 for i in range(len(x))]42 temp_y_arr = []4344 n = len(yArr)45if n == 0:46 res = x47else:48for item in yArr:49 k = QR.dot(x, item) / QR.dot(item, item)50 temp_y_arr.append(QR.list_multi(-k, item))51 temp_y_arr.append(x)5253for item in temp_y_arr:54for i in range(len(item)):55 res[i] += item[i]56return res5758 @staticmethod59def normal(matrix):60 yArr = []61 yArr.append(matrix[0])6263for i in range(1, len(matrix)):64 yArr.append(QR.one_item(matrix[i], yArr))65return yArr6667 @staticmethod68def normalized(lt):69 res = []70 sm = 071for item in lt:72 sm += math.pow(item, 2)73 sm = math.sqrt(sm)74for item in lt:75 res.append(item / sm)76return res7778 @staticmethod79def matrix_T(matrix):80 mat = copy.deepcopy(matrix)81 m = len(mat[0])82 n = len(mat)83for i in range(m):84for j in range(n):85if i < j:86 temp = mat[i][j]87 mat[i][j] = mat[j][i]88 mat[j][i] = temp89return mat9091 @staticmethod92def matrix_multi(mat1, mat2):93 res = []94 rows = len(mat1[0])95 cols = len(mat1)96for i in range(rows):97 temp = [0 for i in range(cols)]98 res.append(temp)99100for i in range(rows):101for j in range(cols):102 sm = 0103for k in range(cols):104 sm += (mat1[k][i] * mat2[j][k])105 res[j][i] = sm106return res107108def execute(self):109 xArr = []110for i in range(self.degree):111 xArr.append(self.get_col(i))112 yArr = QR.normal(xArr)113 self.Q = []114for item in yArr:115 self.Q.append(QR.normalized(item))116117 self.R = QR.matrix_multi(QR.matrix_T(self.Q), xArr) 118return (self.Q, self.R)119120121# A = [122# [1, 0, -1, 2, 1],123# [3, 2, -3, 5, -3],124# [2, 2, 1, 4, -2],125# [0, 4, 3, 3, 1],126# [1, 0, 8, -11, 4]127# ]128# A = [129# [1, 2, 2],130# [2, 1, 2],131# [2, 2, 1]132# ]133 A = [134 [3, 2, 4],135 [2, 0, 2],136 [4, 2, 3]137 ]138139 temp = copy.deepcopy(A)140 val = [] # 特征值141 times = 20 # 迭代次数142for i in range(times):143 qr = QR(temp)144 (q, r) = qr.execute()145 temp = QR.matrix_multi(r, q)146 temp = QR.matrix_T(temp)147148for i in range(len(temp)):149for j in range(len(temp[0])):150if i == j:151 val.append(temp[i][j])152# 特征值153print(val)结果展⽰总结使⽤QR分解迭代求特征值,收敛的⽐较快,也可以求出所有的特征值,但是如果要求特征向量的话,还是需要求解线性⽅程组(感觉很⿇烦)。
求解矩阵特征值问题的算法研究

求解矩阵特征值问题的算法研究求解矩阵特征值问题是线性代数和数值计算中的重要问题。
特征值问题的一般形式为Ax = λx,其中A是一个矩阵,x是一个非零向量,λ是一个标量。
求解特征值问题即寻找矩阵A 的特征值λ和对应的特征向量x。
以下是一些常用的求解矩阵特征值问题的算法:1. 幂迭代法:幂迭代法是求解特征值问题的一种迭代方法。
该方法的基本思想是选择一个随机的初始向量x0,通过迭代计算xk+1 = Axk / ||Axk||,其中||.||表示向量的2-范数。
随着迭代的进行,x收敛到矩阵A的主特征向量,而λ则收敛到对应的主特征值。
2. QR迭代法:QR迭代法是求解特征值问题的一种迭代方法。
该方法首先通过QR分解将矩阵A分解为Q和R的乘积,其中Q是一个正交矩阵,R是一个上三角矩阵。
然后,将分解后的矩阵重新相乘得到新的矩阵A',迭代此过程直到A'的对角线元素收敛到矩阵的特征值。
3. 特征向量分解法:特征向量分解法是求解特征值问题的一种直接方法。
该方法通过对矩阵A 进行特征向量分解得到矩阵V和对角矩阵Λ,其中V的列向量是A的特征向量,Λ的对角线元素是A的特征值。
特征向量分解法在理论上可以得到A的所有特征值和特征向量,但实际计算中对大型矩阵会较为困难。
4. Davidson方法:Davidson方法是求解特征值问题的一种迭代方法。
该方法采用子空间迭代的方式逐步构建特征子空间,并寻找特征值和对应的特征向量。
Davidson方法可以有效地处理大型矩阵,特别适合于求解特征值问题中的一部分特征对。
这些算法都有各自的优点和适用范围,研究者可根据具体问题选择合适的算法进行求解。
矩阵分析与计算QR算法

H11 0
H12 ⎟⎞ p H22 ⎠ n − p
p n− p
这样就可以将问题分解为较小型的问题.如果这种
情况出现在p =n−1或n−2情形,原特征值问题就可以
进行收缩.在实际计算中要对次对角元进行判断,
某个次对角元适当小时就进行分解或收缩.
在定理4.3的条件(4.7)′下,对H的QR迭代中,收敛
(3) Ak = Qk Rk , 其中 Rk = Rk " R2R1
证明 易证(1).从它递推得
Ak+1 = QkT AkQk = (Q1Q2 "Qk )T A(Q1Q2 "Qk )
即得(2),并知Ak+1与A相似,再由 Qk Rk = Q1"Qk Rk " R2 R1 = Qk−1 Ak Rk−1
⎜⎛ 0.6 0.8 0⎟⎞ G(1,2,θ1) = ⎜ − 0.8 0.6 0⎟
⎜⎝ 0 0 1⎟⎠
⎜⎛ 1
G(2,3,θ2 ) = ⎜ 0
最后得到
⎜⎝ 0
0 0.9996 − 0.0249
0 ⎟⎞ 0.0249 ⎟ 0.9996 ⎟⎠
⎜⎛ 4.7600 H2 = ⎜ 0.3200
⎜⎝ 0.0000
基本收敛的概念并未指出{Ak}严格上三角部分 元素是否收敛.但对求A的特征值而言,基本收 敛已足够了.
定理4.3 设矩阵 A∈Rn×n ,其特征值满足
|λ1| > |λ2| > … > |λn| >0 (4.7)′ λi对应特征向量xi ,i =1, 2,…,n.以xi为列的方阵 记为X = (x1, x2 , … , xn ).设X −1可分解为X −1=LU, 其中L为单位下三角阵,U为上三角阵.则QR算法
qr迭代法求特征值特征向量

qr迭代法求特征值特征向量QR迭代法是一种求解矩阵特征值和特征向量的迭代算法。
它基于矩阵的QR分解,通过迭代得到一个上三角矩阵,该矩阵的对角线上的元素就是原矩阵的特征值,而对应的特征向量可以通过反向迭代求得。
QR分解是将一个矩阵分解为一个正交矩阵和一个上三角矩阵的乘积的过程。
对于一个实对称矩阵A,可以通过QR分解得到A=QR,其中Q是一个正交矩阵,R是一个上三角矩阵。
在QR迭代法中,我们通过反复进行QR分解,得到一系列上三角矩阵,最终收敛到一个上三角矩阵,其对角线上的元素就是原矩阵的特征值。
QR迭代法的具体步骤如下:1.将原矩阵A进行QR分解,得到A=QR。
2.计算RQ,得到新的矩阵A'。
3.将矩阵A'再进行QR分解,得到A'=Q'R'。
4.重复步骤2和3,直到矩阵A'收敛到一个上三角矩阵。
5.上三角矩阵的对角线上的元素就是原矩阵的特征值。
在QR迭代法中,每一次迭代的计算量较大,因为需要进行两次矩阵乘法。
为了降低计算量,可以使用对角线P的元素来近似表示特征向量。
具体方法是对矩阵A应用反向迭代,即计算A*x=b,其中x是特征向量的近似表示,b是一个长度为n的列向量,其元素全部为0,只有第i个元素为1。
通过迭代计算,可以得到逼近特征向量的值。
QR迭代法的优点是可以求解任意矩阵的特征值和特征向量,且收敛速度较快。
它不需要知道矩阵的特征值的个数和大小,适用于实对称矩阵和复矩阵。
但它对于存在特征值重复的矩阵求解时收敛速度会较慢,需要进行更多的迭代运算。
最后,我想提醒一下,QR迭代法是一种近似算法,无法保证得到的特征值和特征向量是精确解。
在实际应用中,如果需要得到精确解,可以使用其他更精确的算法,如雅可比迭代法或Hessenberg矩阵迭代法。
普通矩阵特征值的QR算法

普通矩阵特征值的QR 算法摘 要求矩阵的特征值有多种不同的办法,本文主要介绍用QR 法求矩阵的特征值,QR 法是目前求中等大小矩阵全部特征值的最有效方法之一,使用于求实矩阵或复矩阵的特征值,它和雅可比法类似,也是一种变换迭代法。
关键词:QR 分解 迭代序列 特征值 Matlab一 、QR 方法的理论:对任意一个非奇异矩阵(可逆矩阵)A ,可以把它分解成一个正交阵Q 和一个上三角阵的乘积,称为对矩阵A 的QR 分解,即A=QR 。
如果规定R 的对角元取正实数,这种分解是唯一的。
若A 是奇异的,则A 有零特征值。
任取一个不等于A 的特征值的实数μ,则A-μI 是非奇异的。
只要求出A-μI 的特征值和特征向量就容易求出矩阵A 的特征值和特征向量,所以假设A 是非奇异的,不是一般性。
设A=A 1 ,对A 1 作QR 分解,得A 1 = Q 1R 1,,交换该乘积的次序,得A 2 = R 1Q 1=,由于Q 1正交矩阵,A 1到A 2的变换为正交相似变换,于是A 1和A 2就有相同的特征值。
一般的令A 1=A ,对k=1,2,3,…..⎩⎨⎧==+)()(1迭代定义分解kk k k k k Q R A QR R Q A这样,可得到一个迭代序列{A k },这就是QR 方法的基本过程。
二、QR 方法的实际计算步骤Householder A Hessenberg B -----→=用阵作正交相似变换上第阵一步............*::::***⎛⎫⎪ ⎪ ⎪* ⎪** ⎪ ⎪ ⎪ ⎪*⎝⎭Householder 变换:如果 v 给出为单位向量而 I 是单位矩阵,则描述上述线性变换的是 豪斯霍尔德矩阵 (v * 表示向量 v 的共轭转置)H=I -2VV*1k k kGiven k k k B Q R B B R Q +=⎧-------→→⎨=⎩用变换产生迭代序列第二步12***n λλλ⎡⎤⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦Householder A B -----→=用阵作正交相似变换(对称阵)三对角阵*****⎛⎫ ⎪ ⎪⎪*⎪**⎪ ⎪* ⎪ ⎪*⎝⎭三、化一般矩阵为Hessenberg 阵称形如11121112122212323331n n n n n nn nn h h h h h h h h h h h H h h ---⎡⎤⎢⎥⎢⎥⎢⎥=⎢⎥⎢⎥⎢⎥⎣⎦的矩阵为上海森堡(Hessenberg )阵。
QR基本法和位移QR法矩阵特征值求解

实验一:编程实现以下科学计算算法,并举一例应用之。
QR基本法和位移QR法矩阵特征值求解1.QR基本法算法说明:QR基本算法是求矩阵特征值的最有效和应用最广泛的一种方法方法,其基本依据是以下两个定理:1)设A是n阶矩阵,其n个特征值为λ1、λ2、…λm,那么存在一个酉矩阵U使得U T AU是以λ1、λ2、…λm为对角元的上三角矩阵。
2)设A是n阶实矩阵,那么,存在一个正交矩阵Q,使得Q T AQ为一个准上三角矩阵,它的每一个对角元是A的一个特征值,对角元上的二阶块矩阵的两个特征值是A的一对共轭复特征值。
QR基本算法的过程如下:给定循环步数M,A T=A,k=1,2,…M,计算:A k=Q k R kA k+1=R k Q kQR基本算法有如下的收敛性质:如果A的特征值满足|λ1|>|λ2|>|λ2|≥…≥|λm|,则QR基本算法产生的矩阵序列{A k}基本收敛到上三角矩阵(特别,当A为对称阵时,收敛到对角阵),对角元素收敛到A的特征值。
在MATLAB中变成实现的QR基本算法的函数为:qrtz功能:QR基本算法求矩阵全部特征值。
调用格式:l=qrtz(A,M).其中,A为已知矩阵;M为迭代步数;L为矩阵A的全部特征值。
QR基本算法的流程图:l=qrtz(A,M)QR基本算法的MATLAB程序代码如下:function l=qrtz(A,M)for i=1:M[q,r]=qr(A);A=r*q;l=diag(A);endtask11.mformat longA=[1,5,6;4,7,0;8,11,3]l=qrtz(A,20)disp('¾«È·½â')l=eig(A)运行过程和结果:2.位移QR算法位移QR法是为了加快QR算法的收敛速度,其算法的迭代过程如下:给定循环步数M,A1=Hessenberg(A),k=1,2,…M,选择μk,然后计算:A k-μk I=Q k R kA k+1= Q k Q k+μk I一般μk的选择有以下两种考虑方法:(1) 选μk =μk ,即瑞利商位移;(2) 迭代过程中,如果子矩阵的两个特征值为实数时,选最接近a k n,n 的那个作为μk ,即威尔森位移瑞利商位移QR 法流程图如下:在MATLAB 中编程实现的瑞利商位移的QR 算法的函数为:rqrtz 。
- 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。
1
问题描述
对于给定的矩阵A(存在n个不相等的实特征值),采用QR方法求其n个特征值。本文中的矩 阵A由文件”gr9 009 00c rg.mm”给出。并以.mm格式描述。
2
2.1
算法描述
得到实对称矩阵A
本文由文件”gr9 009 00c rg.mm”得到了以.mm格式描述的系数矩阵A。A矩阵是900∗900的大 型稀疏对称矩阵。于是,在matlaB中,使用”A=zeros(900,900)”语句生成900∗900的零矩 阵。再按照.mm文件中的描述,分别对第i行、第j列的元素赋对应的值,就生成了系数矩 阵A,并将A存为.mat文件以便之后应用。 由于A矩阵是实对称的900阶矩阵,所以A矩阵有900个实特征值,可以使用QR方法求解。 下面描述具体算法。
2.4
用 吉 文斯 变 换 实 现 一次 等 效 的 QR迭 代
在上文,已用豪斯荷尔德变换将矩阵A化为上海森伯格矩阵AH . 对于一次QR迭代而言,相 当于用吉文斯变换实现等效地实现以下过程: G(n − 1, n; θn−1 )G(n − 2, n − 1; θn−1 ) · · · · · · G(1, 2; θ1 )AH = R 其中G(i − 1, i; θi−1 )是吉文斯扩充矩阵,是正交矩阵。 该过程等效于: AH = QR 其中 Q = G(1, 2; θ1 )T · · · · · · G(n − 2, n − 1; θn−1 )T G(n − 1, n; θn−1 )T 得到上三角矩阵R后,再进行: AHnew = RQ = RG(1, 2; θ1 )T · · · · · · G(n − 2, n − 1; θn−1 )T G(n − 1, n; θn−1 )T
⫼ঢ়᭛ᮃবᤶᅲ ⦄AHⶽ䰉ⱘϔ QR䗁ҷ
Ϟ⍋Ểԃ Ḑⶽ䰉AH ᵘ䗴ঢ়᭛ᮃবᤶG(i,i+1)ᬍবAHⶽ䰉i㸠i+1㸠ⱘ ؐˈՓᕫAH i+1,i =0ˈㄝᬜഄᅲ⦄˖ A=QR.
N O i=1:n-1 YES
T T
ᕫࠄϞϝ㾦 ⶽ䰉R
ᣝ✻G(1,2) ,G(2,3) ĂĂG(n,n-1) ᴹձᬍবRⶽ䰉 ⱘi߫i+1߫ˈㄝᬜഄᅲ⦄˖ Anew=RQDŽ N O i=1:n-1 YES 㒧ᴳᴀᄤᑣˈ ᅲ⦄ձQR䗁 ҷᕫࠄAnew
3
结果输出
在f unction[eigenvalue, f lag − 1] = QR(A, n)函数中,设置accuracy = 10−3 ,取n=900,A为 本文2.1所得的矩阵,运行该函数。结果显示,f lag1 = 9344。也就是,经过9344次迭代, 上海森伯格矩阵次AH 对角元素接近0,可由对角元素得到A的特征值。所得结果在同一文件 夹中的eigenvalue.mat中 图四显示了900个特征值从小到大的排列。 表一给出了900个特征值中最大的八个和最小的八个特征值。 6
QR迭代法求解矩阵A的特征值
沈欢00986096 北京大学工学院,北京100871 2011 年 10 月 23 日 摘 要
本文用QR迭代法求解矩阵A的特征值:第一步先用豪斯荷尔德变换将矩阵A化为上海 森伯格矩阵AH ,第二步再对AH 进行QR迭代(使用吉文斯变换),当迭代满足精度要求时 输出n个特征值。
Yes 䖒ࠄQR䗁ҷⱘ㊒ᑺ˛
No
ᕫࠄ ⒵䎇 ㊒ᑺ ⱘϞ ϝ㾦 ⶽ䰉 An
㒧 ᴳ ᑣ
⫼ঢ়᭛ᮃবᤶᅲ⦄ϔ;1=RiQi
图 1: 主程序(QR迭代法求矩阵特征值)算法 用豪斯荷尔德变换将A矩阵化为上海森伯格矩阵AH 是实现下述过程: A1 = H1 AH1 A2 = H2 A1 H2 ······ An−2 = Hn−2 An−3 Hn−2 其中H1 、H2 · · · · · · 满足: H1 a11 a21 a31 . . . an1 a11 a21 0 . . . 0
=
(2)
3
H2
a12 a22 a32 a42 . . . an2
=
a12 a22 a32 0 . . . 0
(3)
······ 于是构造豪斯荷尔德矩阵Hi 成为此函数的核心。对于Hi ,目的是将之前得到的Ai−1 中的 第i列实现如下变换: a1i a1i a2i a2i a3i a3i ··· ··· a a (4) Hi ii ii = ai+1,i a i+1,i ai+2,i 0 . . . . . . 0 ani 于是从Ai−1 中取出n − i维列向量 − → x = [ai+1,i , ai+2,i , · · · , ani ]T 取 → α= − x 2; 和n-i维单位向量 构造 − → e = [1, 0, 0, · · · , 0]T − → → x − α− e ω= − → − → x −αe 2
2.2
QR算 法
QR算法的核心思想是:对一个具有n个不相等的特征值的矩阵A,进行如下变换: A = Q1 R1 ; A1 = R1 Q1 A1 = Q2 R2 ; A2 = R2 Q2 ········· 1
Ak = Qk+1 Rk+1 ; Ak+1 = Rk+1 Qk+1 ········· 其中Qi 是正交矩阵,Ri 是可逆的上三角矩阵。由于:
于是就得到了(n − i) ∗ (n − i)的豪斯荷尔德矩阵Hi , 再将其扩充: Hi = 得到Hi 。实现: Ai = Hi Ai−1 Hi 将此步依次从i = 1至i = n − 2操作n-2次,得到上海森伯格矩阵AH ,且与原矩阵A相似。 本函数(function AH=A-to-AHessenberg(A,n))流程图如图二所示。 4 Ii∗i 0 0 Hi (5)
2.3
用 豪 斯 荷尔 德 变 换 将A矩 阵化 为 上 海 森伯 格 矩 阵 AH
H = I − 2wwT
豪斯荷尔德矩阵 其中w是单位向量。所以H矩阵是正交、对称矩阵。
2
⫼QR䗁ҷ⊩∖ Aⱘ⡍ᕕؐ ˄Џᑣ˅
A ⶽ 䰉
⫼䈾ᮃ㥋ᇨᖋবᤶᇚAবЎϞ⍋ỂԃḐⶽ 䰉A H 䕧ߎ An ⱘᇍ 㾦ܗ ㋴ ЎA ⶽ䰉 ⱘ⡍ ᕕؐ
2 AHi,i 2 + AHi +1,i
si = 于是,第i行和第i+1行的元素变换为:
AHi+1,i
2 + AH 2 AHi,i i+1,i
AHi,j = ci ∗ AHi,j + si ∗ AHi+1,j AHi+1,j = ci ∗ AHi+1,j − si ∗ AHi,j j = 1, 2, · · · · · · , n . 对i从1至n-1,重复以上过程,则等效实现QR分解,并得到上三角矩阵R。 在第二个过程中,R矩阵右乘上Q矩阵,即: AHnew = RQ = RG(1, 2; θ1 )T · · · · · · G(n − 2, n − 1; θn−1 )T G(n − 1, n; θn−1 )T 等效于依次对R矩阵的第i列和第i+1列进行列变换。变换后的元素为: Rj,i = ci ∗ Rj,i + si ∗ Rj,i+1 Rj,i+1 = ci ∗ Rj,i+1 − si ∗ Rj,i j = 1, 2, · · · · · · , n . 对i从1至n-1,重复以上过程,则等效地实现AHnew = RQ过程。 至此,完成了一次QR迭代。 本函数(AHnew=QR-iteration-once(AH,n))的流程图如图三所示。
⫼䈾ᮃ㥋ᇨ ᖋবᤶᇚAⶽ 䰉ⳌԐবᤶ ЎϞ⍋Ểԃ Ḑⶽ䰉AH i=n-2?
A ⶽ 䰉 䕧ߎ Ϟ⍋ Ểԃ Ḑⶽ 䰉AH
Ai-1 ᵘ䗴n-i㓈ⱘwi䞣
Ai
߽⫼wi ᵘ䗴(n-i)*(n-i)ⱘH(n-i)*(n-i)
ᠽܙH(n-i)*(n-i)ᕫࠄHi
Ai=Hi*Ai-1*Hi
图 2: 用豪斯荷尔德变换将矩阵A化为上海森伯格矩阵AH算法
T Ai+1 = QT i+1 (Qi+1 Ri+1 )Qi+1 = Qi+1 Ai Qi+1
所以Ai+1 和Ai 相似,具有相同的特征根。当i→ ∞时,Ai 对角线以下元素趋于零。 利用这一原理,可以对给定A进行QR迭代,最终得到上三角矩阵(对角元素不一定非零), 输出对角元素可以作为A的特征值。 我们又发现,具有如下上海森伯格形式的矩阵AH : a11 a12 a13 . . . a1n−1 a1n a21 a22 a23 . . . a2n−1 a2n 0 a32 a33 a34 ... ... (1) AH = ... ... ... ... 0 0 . . . . . . an,n−1 ann AH 在QR迭代过程中可以保持上海森伯格形式不变。 所以,QR迭代算法分两步实现,第一步将矩阵A变换为上海森伯格矩阵AH ,第二步对此上 海森伯格矩阵AH 进行QR迭代,直至满足精度为止,输出结果。 在第一步中,将矩阵A变换为上海森伯格矩阵AH ,为了保持变换后的AH 与A具有相同的特 征值,这里采用豪斯荷尔德变换(正交变换),本步在后文中详细说明;在第二步中,对上 海森伯格矩阵AH 进行QR迭代,并不去直接找出Q矩阵和R矩阵,而是采用吉文斯变换达到 相同的效果,本步在后文中也有详细说明;对于QR迭代的精度,在算法采用:当迭代到 第k次时有: |Ai+1,i | < accuracy |Ai,i | + |Ai+1,i+1 | 对于i=1,2,· · · · · · n-1成立,其中accuracy是所设精度。认为此时的AH 次对角元素已经为 零,迭代结束,输出该矩阵对角元素作为A的特征值。 QR算法主程序(见程序:function [eigenvalue,flag1]=QR(A,n))的程序流程图如图1所示。 对于方框中的两步:1、用豪斯荷尔德变换将A矩阵化为上海森伯格矩阵AH ; 2、用吉文斯 变换实现一次等效的QR迭代;将分别描述。