《計算方法》課件第8章

第8章矩陣特征值及特征向量的數(shù)值求解8.1引言8.2冪法與反冪法8.3雅可比方法8.4

QR方法

8.1引言

求矩陣的特征值和特征向量是代數(shù)計算中的重要課題，同時在自然科學和工程計算中，很多問題也常常歸結(jié)為求解矩陣的特征值和特征向量的問題。

根據(jù)線性代數(shù)理論課程所講述的計算矩陣特征值和相應特征向量的方法，若已知n階矩陣A：求矩陣A的特征值就是求A的特征方程的解，其特征方程為即(8.1)式(8.1)是以λ為未知數(shù)的一元n次方程，通過解特征方程可得到矩陣A的特征值λ1，λ2，…，λn

。關于矩陣的特征向量，則通過如下的方程進行求解：(8.2)式(8.2)實質(zhì)為方程組，通過計算可得到λi對應的特征向量xi(i=0，1，…，n)。當矩陣的階數(shù)n較小時，可以利用上述方法求解矩陣的特征值和特征向量，而當n較大時，則存在困難。另外,實際問題中對求解矩陣的特征值和特征向量有不同的具體要求，如要求計算按模(絕對值)最大的特征值和特征向量或全部特征值和特征向量等。為此，本章介紹常用的求解矩陣特征值和特征向量的數(shù)值方法。 8.2冪法與反冪法

8.2.1冪法

1.冪法的基本原理

冪法是一種求實矩陣按模最大的特征值和相應特征向量的迭代方法。它是通過迭代產(chǎn)生向量序列，由此計算矩陣的特征值和特征向量。

設n階實方陣A有n個線性無關的特征向量x1，x2，…，xn，其所對應的特征值為λ1，λ2，…，λn，并且已按模的大小排列，即(8.3)冪法的目標是求解方陣A的特征值λ1及其特征向量x1，下面介紹其基本原理。

首先任意給定非零初始向量，設為v0≠0，并給出迭代公式：(8.4)依據(jù)該迭代公式可以得到一個向量序列V={v0，v1，…，vn}。由于x1，x2，…，xn線性無關，構(gòu)成Rn上的一組基，因此對于任意給定的非零向量v0必然存在一組不全為0的數(shù)a1，a2，…，an，使得將v0代入式(8.4)有(8.5)根據(jù)線性代數(shù)的知識有所以根據(jù)式(8.5)有(8.6)同理可得(8.7)

迭代一直進行下去，可得到一個向量序列V={v0，v1，…，vn}。對其迭代過程中的向量進行分析可知，當?shù)降趉步和第k+1步時有（8.9）由于 ，則因此,當k足夠大時，根據(jù)式(8.8)和式(8.9)有進而得到(8.10)式(8.10)說明當k足夠大時，向量vk+1和vk近似地為線性關系，常數(shù)λ1即是實方陣A的按模最大的特征值，vk+1和vk

是對應的特征向量，具體求λ1時，可取(8.11)其中(vk+1)j和(vk)j分別表示vk+1和vk的第j個分量。上述內(nèi)容即為冪法計算按模最大的特征值和相應特征向量的過程，可以看出,該過程主要是計算矩陣A的冪Ak及其與已知向量v0的乘積Akv0，因此被稱為乘冪法，簡稱冪法。其實質(zhì)是一種迭代方法，迭代的收斂速度主要取決于比值|λ2/λ1|，|λ2/λ1|越小，收斂速度越快。

2.冪法的歸一化處理

在冪法的計算過程中，每一步迭代之前都進行歸一化處理，由此得到歸一化后的迭代公式：(8.12)由式(8.12)可得出如下兩個結(jié)論：(2)

證明

(1)根據(jù)式(8.12)有而所以故由于故結(jié)論（1）得證。

(2)根據(jù)式(8.12)有所以結(jié)論(2)得證。

例8－1

已知矩陣用歸一化后的冪法迭代公式計算A按模最大的特征值及相應的特征向量，要求

解取v0=(1，1，1)T，根據(jù)迭代公式可得…表8.1例8－1的計算結(jié)果

當?shù)羕=9時，|max(v9)－max(v8)|=0.0004147221<0.001，達到了精度要求，停止迭代。由此可得矩陣A的按模最大特征值λ1近似為44.999931，對應的特征向量為u9=［0.000131，0.499868，1.000000］T。8.2.2反冪法

1.反冪法的基本原理

根據(jù)冪法求矩陣按模最大的特征值及其相應特征向量的思想可得到求矩陣按模最小的特征值及其相應特征向量的方法，即反冪法，其基本原理如下。

設A為n階非奇異矩陣，有n個線性無關的特征向量x1，x2，…，xn，其對應的特征值為λ1，λ2，…，λn，且特征值滿足

2.反冪法的計算公式

根據(jù)反冪法的基本原理，結(jié)合冪法的計算公式(8.12)可得到反冪法的計算公式如下：（8.14）與式(8.12)相同，根據(jù)式(8.14)也可得到如下的兩個結(jié)論：（1）（2）根據(jù)式(8.14)可通過迭代計算得到矩陣A的按模最小的特征值及其相應的特征向量。從式(8.14)可以看出，必須計算A矩陣的逆矩陣A－1，而計算逆矩陣往往比較困難，因此將式(8.14)中的vk=A－1uk－1改為Avk=uk－1，則式(8.14)又可寫為（8.15）（8.16）

例8－2用反冪法計算例8－1中矩陣的按模最小的特征值及相應的特征向量。

解對矩陣A做LU分解，有于是，根據(jù)式(8.16)取初始向量u0=v0=［1，1，1］T，計算過程略，計算結(jié)果如表8.2所示。表8.2例8－2的計算結(jié)果

8.3雅可比方法

雅可比方法是求實對稱矩陣全部特征值及對應特征向量的方法。它也是一種迭代法，其基本思想是把對稱矩陣A經(jīng)一系列正交相似變換化為一個近似對角矩陣，該對角矩陣的對角元素就是A的近似特征值，對應的特征向量可由各個正交變換矩陣的乘積求得。

雅可比方法是求實對稱矩陣全部特征值及對應的特征向量的方法。它也是一種迭代法，其基本思想是把對稱矩陣A經(jīng)一系列正交相似變換約化為一個近似對角陣，從而該對角陣的對角元就是A的近似特征值，由各個正交變換陣的乘積可得對應的特征向量。

先從二階矩陣開始討論，設二階實對稱陣為尋求正交矩陣使得A經(jīng)過正交相似變換化為對角矩陣。設正交矩陣R為則(8.17)其中,因為R是正交矩陣，故B是A的正交相似矩陣，為使B成為對角矩陣，只須適當選取θ，使即滿足：(8.18)由此可以確定θ，從而確定正交矩陣R，由此得A的特征值為(8.19)由，可以得：如果記R的列向量為和則有：其中λi為A的特征值，xi就是對應的特征向量，根據(jù)R可得A的兩個特征向量為(8.20)將二階矩陣推廣到n階矩陣的情況，設矩陣A∈Rn×n是對稱矩陣，對A作正交旋轉(zhuǎn)相似變換，取n階正交旋轉(zhuǎn)矩陣R為ij(8.21)矩陣A經(jīng)過正交相似變換后得到矩陣A1，將矩陣A1與矩陣A比較可得：A1的第i行、第j行、第i列及第j列的元素發(fā)生了變化，其他元素不變，即(8.22)如果取θ使得，即滿足(8.23)可得

定理8.1當k→∞時，Ak→Λ=diag(λ1，λ2，…，λn)，其中λ1，λ2，…，λn即為矩陣A的特征值。

證明根據(jù)式(8.22)將矩陣A1和A比較可得

(1)(8.24)即經(jīng)過正交旋轉(zhuǎn)變換后，矩陣所有元素的平方和不變。(2)(8.25)利用式(8.23)，式(8.25)可變?yōu)?8.26)根據(jù)式(8.24)和式(8.26)可知，矩陣A經(jīng)過正交旋轉(zhuǎn)相似變換后，矩陣中的對角線元素的平方和比原來增加了2a2ij，而非對角線元素的平方和減少了2a2ij

。如果取則有（8.27）其中，S為A中非對角線元素的平方和。則矩陣A通過一次變換后得到矩陣A1，A1的非對角線元素的平方和S1為(8.28)將式(8.27)代入式(8.28)可得(8.29)當k

時，

即矩陣非對角線元素的平方和趨于0，也就是非對角線元素趨于0，所以(8.30)令

例8－3用雅可比方法求對稱矩陣的特征值及特征向量，要求

解第一步:找出矩陣A中非對角元素中絕對值最大的元素，取a12=－1，則i=1，j=2由式(8.23)可得由式(8.23)可得則所以可得正交旋轉(zhuǎn)矩陣：令A1=A，得一次正交相似變換后的A2為RT1A1R1，即同時有可以看出,A2中非對角線元素的絕對值的最大值不滿足精度要求。第二步:與第一步的計算過程相同。

直到第n步，滿足所要求的精度。從第二步開始，計算過程略，計算結(jié)果如表8.3所示。當?shù)恋?步時，A6中非對角線元素的絕對值的最大值為0.016758，此時0.016758<0.05，達到精度要求，迭代停止?！捎贏6與A有相同的特征值，故A的特征值近似為對應的特征向量分別為表8.3例8-3的計算結(jié)果 8.4

QR方法

8.4.1

QR方法的基本思想

QR方法是求解矩陣全部特征值的一種有效方法，也是一種變換迭代法。下面介紹QR方法的基本思想。

設n階矩陣A為非奇異矩陣，通過一系列正交變換，得到A的正交三角分解(8.31)其中，Q是正交矩陣，R是上三角形矩陣。設A1=A，對A1進行QR分解有(8.32)若令A2=R1Q1，由式(8.32)可得故可得(8.33)由于Q1是正交矩陣，故有QT1=Q－11，進而有(8.34)則A1與A2相似，具有相同的特征值。對比A1與A2發(fā)現(xiàn),A2經(jīng)過A1正交變換，相對A1更接近于一個上三角矩陣，且求矩陣的特征值更容易一些，進而可將求A1的特征值轉(zhuǎn)化為求A2的特征值。繼續(xù)對A2進行QR分解有A2=Q2R2，令A3=R2Q2，從而有A3=QT2A2Q2，求A2的特征值轉(zhuǎn)化為求A3的特征值。如此繼續(xù)下去，到第k步時得到Ak，對Ak進行QR分解有(8.35)并令Ak+1=RkQk，則(8.36)則求Ak特征值轉(zhuǎn)化為求Ak+1的特征值。從A1到Ak+1，逐漸接近于一個上三角矩陣，由于A1,A2,A3,…,Ak,Ak+1有相同的特征值，若迭代的次數(shù)足夠多，則Ak+1近似為一個上三角矩陣，那么求原矩陣A的特征值就轉(zhuǎn)化為求一個上三角矩陣特征值的問題。8.4.2矩陣的QR分解

設A為n階矩陣，下面介紹矩陣的QR分解過程。

設平面旋轉(zhuǎn)矩陣R(i，j)為ij將R(i，j)左乘A，可將A的第i列第j行的元素化為零。記R(i，j)A=［a(1)ij］n×n，有(8.37)對于式(8.37)，如果θ滿足，則有a(1)