4.2 矩陣特征值和奇異值
對于n階方陣A，求數(shù)λ和向量X，使得等式AX=λX成立，滿足等式的數(shù)λ稱為A的特征值，向量X稱為A的特征向量。方程AX=λX和(A-λI)X=0是兩個(gè)等價(jià)方程，要使方程(A-λI)X=0有非0解X，則必須使其行列式等于0，即|A-λI |=0。
由線性代數(shù)可知，行列式|A-λI |是一個(gè)關(guān)于λ的n階多項(xiàng)式，因此方程|A-λI |=0是一個(gè)n次方程，有n個(gè)根（包含重根）。n個(gè)根就是矩陣A的n個(gè)特征值，每一個(gè)特征值對應(yīng)無窮多個(gè)特征向量。所以，矩陣的特征值問題有確定的解，但特征向量問題沒有確定的解。
4.2.1 特征值和特征向量的求取
特征值和特征向量在科學(xué)研究和工程計(jì)算中的應(yīng)用非常廣泛。在MATLAB中，計(jì)算矩陣A的特征值和特征向量的函數(shù)是eig(A)，常用的調(diào)用格式有以下3種。
(1)E=eig(A)：用于求矩陣A的全部特征值，構(gòu)成向量E。
(2)[V,D]=eig(A)：用于求矩陣A的全部特征值，構(gòu)成對角矩陣D，并求A的特征向量，構(gòu)成V的列向量。
(3)[V,D]=eig(A,’nobalance’)：與上一種格式類似，只是上一種格式中是先對A作相似變換后再求矩陣A的特征值和特征向量，而本格式中則是直接求矩陣A的特征值和特征向量。
一個(gè)矩陣A的特征向量有無窮多個(gè)，eig函數(shù)只找出其中的n個(gè)，A的其他特征向量均可由這n個(gè)特征向量組合表示。
【例4-11】 簡單實(shí)數(shù)矩陣的特征值示例。
>> A=[1,-3;2,2/3]
A =
1.0000 -3.0000
2.0000 0.6667
>> [V,D]=eig(A)
V =
0.7746 + 0.0000i 0.7746 + 0.0000i
0.0430 - 0.6310i 0.0430 + 0.6310i
D =
0.8333 + 2.4438i 0.0000 + 0.0000i
0.0000 + 0.0000i 0.8333 - 2.4438i
【例4-12】 矩陣中有元素與截?cái)嗾`差相當(dāng)時(shí)的特性值問題。
>> A=[3 -2 -0.9 2*eps
-2 4 -1 -eps
-eps/4 eps/2 -1 0
-0.5 -0.5 0.1 1 ];
>> [V1,D1]=eig(A); % 求特征值
>> ER1=A*V1-V1*D1 % 查看計(jì)算誤差
ER1 =
0.0000 0 -0.0000 0.0000
0 -0.0000 0.0000 -0.0000
0.0000 -0.0000 -0.0000 0
0 0.0000 0.0000 -0.3789
>> [V2,D2]=eig(A,'nobalance');
>> ER2=A*V2-V2*D2 % 查看計(jì)算誤差
ER2 =
1.0e-14 *
-0.1776 -0.0111 -0.0559 -0.0833
0.3553 0.1055 0.0343 0.0555
0.0017 0.0002 0.0007 0
0.0264 -0.0222 0.0222 0.0333
在本例，若求特征值的過程中不采用'nobalance'參數(shù)，那么計(jì)算結(jié)果是具有相當(dāng)大的計(jì)算誤差的。這是因?yàn)樵趫?zhí)行eig命令的過程中，首先要調(diào)用使原矩陣各元素大致相當(dāng)?shù)摹捌胶狻背绦?，這些“平衡”程序使得原來方陣中本可以忽略的小元素（本例中如eps）的作用被放大了，所以產(chǎn)生了較大的計(jì)算誤差。
但是這種誤差被放大的情況只發(fā)生在矩陣中有元素與截?cái)嗾`差相當(dāng)?shù)臅r(shí)候，在一般情況下，“平衡”程序的作用是減小計(jì)算誤差。
【例4-13】 函數(shù)eig 與eigs 的比較示例。eigs函數(shù)是計(jì)算矩陣最大特征值和特征向量的函數(shù)。
>> rand('state',0) % 設(shè)置隨機(jī)種子，便于讀者驗(yàn)證
>> A=rand(100,100)-0.5;
>> t0=clock;[V,D]=eig(A);T_full=etime(clock,t0) % 函數(shù)eig的運(yùn)行時(shí)間
>> options.tol=1e-8; % 為eigs設(shè)置計(jì)算精度
>> options.disp=0; % 不顯示中間迭代結(jié)果
>> t0=clock;[v,d]=eigs(A,1,'lr',options); % 計(jì)算最大實(shí)部特征值和特征向量
>> T_part=etime(clock,t0) % 函數(shù)eigs的運(yùn)行時(shí)間
>> [Dmr,k]=max(real(diag(D))); % 在eig求得的全部特征值中找實(shí)部最大的
>> d,D(1,1)
運(yùn)行結(jié)果為：
T_full =
0.9510
T_part =
1.6540
d =
2.7278 + 0.3006i
ans =
2.5933 + 1.5643i
>> vk1=V(:,k); % 與d相同的特征向量應(yīng)是V的第k列
>> vk1=vk1/norm(vk1);v=v/norm(v); % 向量長度歸一
>> V_err=acos(norm(vk1'*v))*180/pi % 求復(fù)數(shù)向量之間的夾角
V_err =
8.5377e-007
>> D_err=abs(D(k,k)-d)/abs(d) % 求兩個(gè)特征值間的相對誤差
D_err =
2.6420e-009
在本例中，對函數(shù)運(yùn)行所需要的時(shí)間進(jìn)行了評估。需要指出的是：在實(shí)際使用中因?yàn)橛?jì)算機(jī)的配置和系統(tǒng)狀態(tài)不同，評估得到的絕對時(shí)間也不盡相同，不過我們可以通過同一臺計(jì)算機(jī)上兩種函數(shù)運(yùn)行所需要時(shí)間比較得到兩種算法的優(yōu)劣。通過本例可以得出結(jié)論：使用eigs函數(shù)求一個(gè)特征值和特征向量所需要的時(shí)間，反而比使用eig函數(shù)求全部特征值和特征向量的時(shí)間多。
【例4-14】 用求特征值的方法，求解方程的根。
求解方程組要先構(gòu)造與方程對應(yīng)的多項(xiàng)式的伴隨矩陣A，再求A的特征值，伴隨矩陣A的特征值即為方程的解。具體過程如下：
>> B=[1,-2,3,4,-5]
B =
1 -2 3 4 -5
>> A=compan(B) % 求B 的伴隨矩陣
A =
2 -3 -4 5
1 0 0 0
0 1 0 0
0 0 1 0
>> C=eig(A)
C =
1.1641 + 1.8573i
1.1641 - 1.8573i
-1.1973 + 0.0000i
0.8691 + 0.0000i
>> D=roots(B) % 直接求多項(xiàng)式B的零點(diǎn)
D =
1.1641 + 1.8573i
1.1641 - 1.8573i
-1.1973 + 0.0000i
0.8691 + 0.0000i函數(shù)compan計(jì)算的是矩陣B的伴隨矩陣A。伴隨矩陣的特征值C就是方程的根。roots函數(shù)用于直接對線形方程求解，得到結(jié)果D，可以看出兩種方法得出的結(jié)果C和D是一樣的。
4.2.2 奇異值分解
如果存在兩個(gè)矢量u、v及一個(gè)常數(shù)s，使得矩陣A滿足下式：
Av=su
A’u=sv
則稱s為奇異值，稱u、v為奇異矢量。
將奇異值寫成對角方陣∑，而相對應(yīng)的奇異矢量作為列矢量，則可寫成兩個(gè)正交矩陣U、V，使得：
AV=U∑
A’U=V∑
因?yàn)?span>U、V正交，所以可得奇異值的表達(dá)式為：
A=U∑V’
一個(gè)m行n列的矩陣A經(jīng)奇異值分解，可求得m行m列的距陣U，m行n列的矩陣∑，n行n列的矩陣V。
奇異值分解是另一種正交矩陣分解法，是最可靠的分解法，但是它與QR分解法相比，要花近10倍的計(jì)算時(shí)間。奇異值分解由svd函數(shù)實(shí)現(xiàn)，其調(diào)用格式為：[U,S,V]=svd(A)。
【例4-15】 奇異值分解示例。
>> A=magic(4)
A =
16 2 3 13
5 11 10 8
9 7 6 12
4 14 15 1
>> [U,S,V] = svd(A) % 奇異值分解
U =
-0.5000 0.6708 0.5000 -0.2236
-0.5000 -0.2236 -0.5000 -0.6708
-0.5000 0.2236 -0.5000 0.6708
-0.5000 -0.6708 0.5000 0.2236
S =
34.0000 0 0 0
0 17.8885 0 0
0 0 4.4721 0
0 0 0 0.0000
V =
-0.5000 0.5000 0.6708 0.2236
-0.5000 -0.5000 -0.2236 0.6708
-0.5000 -0.5000 0.2236 -0.6708
-0.5000 0.5000 -0.6708 -0.2236
>> U*S*V' % 分解結(jié)果正確性驗(yàn)證
ans =
16.0000 2.0000 3.0000 13.0000
5.0000 11.0000 10.0000 8.0000
9.0000 7.0000 6.0000 12.0000
4.0000 14.0000 15.0000 1.0000
