|
4.2 矩陣特征值和奇異值 對(duì)于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è)特征值對(duì)應(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)成對(duì)角矩陣D,，并求A的特征向量，構(gòu)成V的列向量,。 (3)[V,D]=eig(A,’nobalance’)：與上一種格式類似,，只是上一種格式中是先對(duì)A作相似變換后再求矩陣A的特征值和特征向量，而本格式中則是直接求矩陣A的特征值和特征向量,。 一個(gè)矩陣A的特征向量有無窮多個(gè),，eig函數(shù)只找出其中的n個(gè)，A的其他特征向量均可由這n個(gè)特征向量組合表示,。 【例4-11】 簡(jiǎn)單實(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); % 向量長(zhǎng)度歸一 >> 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è)特征值間的相對(duì)誤差 D_err = 2.6420e-009 在本例中,，對(duì)函數(shù)運(yùn)行所需要的時(shí)間進(jìn)行了評(píng)估。需要指出的是：在實(shí)際使用中因?yàn)橛?jì)算機(jī)的配置和系統(tǒng)狀態(tài)不同，評(píng)估得到的絕對(duì)時(shí)間也不盡相同,，不過我們可以通過同一臺(tái)計(jì)算機(jī)上兩種函數(shù)運(yùn)行所需要時(shí)間比較得到兩種算法的優(yōu)劣,。通過本例可以得出結(jié)論：使用eigs函數(shù)求一個(gè)特征值和特征向量所需要的時(shí)間，反而比使用eig函數(shù)求全部特征值和特征向量的時(shí)間多,。 【例4-14】 用求特征值的方法，求解方程的根,。 求解方程組要先構(gòu)造與方程對(duì)應(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ù)用于直接對(duì)線形方程求解,，得到結(jié)果D,，可以看出兩種方法得出的結(jié)果C和D是一樣的。 4.2.2 奇異值分解 如果存在兩個(gè)矢量u,、v及一個(gè)常數(shù)s,，使得矩陣A滿足下式： Av=su A’u=sv 則稱s為奇異值，稱u,、v為奇異矢量,。 將奇異值寫成對(duì)角方陣∑，而相對(duì)應(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
|
|
|
