簡介這份資源面向?qū)W習(xí)MATLAB數(shù)值計(jì)算與線性代數(shù)的學(xué)生和工程人員聚焦矩陣處理中的特征值、特征向量與對(duì)角化等核心問題適合作為課堂實(shí)驗(yàn)或自學(xué)練習(xí)的配套材料。壓縮包內(nèi)共1個(gè)doc文檔約98KB內(nèi)容以實(shí)驗(yàn)講義形式組織涵蓋特征多項(xiàng)式求解、eig函數(shù)求特征值與特征向量、實(shí)對(duì)稱矩陣正交對(duì)角化、balance函數(shù)提高計(jì)算精度以及jordan函數(shù)求若當(dāng)標(biāo)準(zhǔn)形等知識(shí)點(diǎn)并配有可直接運(yùn)行的MATLAB代碼示例與運(yùn)算結(jié)果對(duì)照。文檔還整理了trace、sum、prod、cdf2rdf等常用矩陣函數(shù)的用法說明便于讀者在練習(xí)中隨時(shí)查閱。目前已有727人學(xué)習(xí)下載可作為矩陣分析入門與實(shí)驗(yàn)報(bào)告的參考幫助讀者理解特征值分解的數(shù)學(xué)原理并掌握相應(yīng)MATLAB命令的調(diào)用方式。1. 矩陣處理到底在算什么從“找方向”說起很多人第一次聽到“特征值與特征向量”腦子里冒出的是一堆符號(hào)和公式完全不知道它在生活里對(duì)應(yīng)什么。其實(shí)可以這樣理解你把一個(gè)矩陣想象成一臺(tái)“變形機(jī)器”任何向量丟進(jìn)去都會(huì)被拉伸、旋轉(zhuǎn)、翻轉(zhuǎn)。但總有那么幾個(gè)特殊方向丟進(jìn)去之后只被拉長或壓短方向紋絲不動(dòng)——這幾個(gè)方向就是特征向量拉伸的倍數(shù)就是特征值。Matlab 矩陣處理里最核心的一塊就是幫你把這幾條“不動(dòng)方向”和“拉伸倍數(shù)”找出來。搞懂這件事的價(jià)值在于矩陣對(duì)角化、若當(dāng)標(biāo)準(zhǔn)形、實(shí)對(duì)稱矩陣的正交分解全都建立在這個(gè)概念之上。你不需要手算行列式Matlab 的eig、jordan、balance這些函數(shù)已經(jīng)把底層運(yùn)算封裝好了。但前提是你得知道每個(gè)函數(shù)在做什么、結(jié)果怎么讀、什么時(shí)候會(huì)翻車。這篇就把這些概念用大白話拆開讓你看完能分清e(cuò)ig和eigs的區(qū)別知道jordan為什么有時(shí)候“不靈”也明白balance到底在平衡什么。2. 特征值與特征向量Matlab 里到底怎么算、結(jié)果怎么讀2.1 從特征多項(xiàng)式到poly和roots在 Matlab 里求特征值最“笨”但也最直觀的路子是先算特征多項(xiàng)式。特征多項(xiàng)式的定義是 det(λI ? A)其中 I 是單位陣λ 是符號(hào)變量。Matlab 提供了poly(A)直接返回多項(xiàng)式系數(shù)的向量形式。比如一個(gè) 2×2 矩陣poly會(huì)給你三個(gè)系數(shù)分別對(duì)應(yīng) λ2、λ1 和常數(shù)項(xiàng)。A [1 0; 2 3]; p poly(A); % 返回特征多項(xiàng)式系數(shù)向量 f poly2str(p, x); % 轉(zhuǎn)成可讀字符串形式poly返回的是數(shù)值向量poly2str只是把它排版成人類能讀的表達(dá)式。如果你想知道特征值具體是多少還得對(duì)多項(xiàng)式求根用roots函數(shù)。這條路適合教學(xué)演示因?yàn)槟隳苡H眼看到特征多項(xiàng)式長什么樣。但實(shí)際工程里沒人這么干——高階矩陣的多項(xiàng)式求根數(shù)值穩(wěn)定性很差Matlab 內(nèi)部也不會(huì)用這條路去算eig。另一種更“數(shù)學(xué)”的寫法是用符號(hào)變量直接算行列式syms x; E eye(2); f det(x*E - A); % 符號(hào)形式的特征多項(xiàng)式這樣得到的是符號(hào)表達(dá)式適合推導(dǎo)和驗(yàn)證但不適合大規(guī)模數(shù)值計(jì)算。符號(hào)運(yùn)算慢而且當(dāng)矩陣階數(shù)上去之后表達(dá)式會(huì)膨脹到?jīng)]法看。2.2eig的三種返回形式與eigs的區(qū)別真正干活用的是eig。它有三種常見調(diào)用方式返回的東西不一樣用途也不一樣。第一種d eig(A)只返回一個(gè)列向量里面裝著所有特征值。你只關(guān)心特征值大小、不關(guān)心特征向量時(shí)用這個(gè)。第二種[V, D] eig(A)返回兩個(gè)矩陣。D 是對(duì)角矩陣主對(duì)角線上的元素就是特征值V 的每一列是對(duì)應(yīng)的特征向量。滿足的關(guān)系是 A·V V·D等價(jià)于 D V?1AV。這就是矩陣對(duì)角化的核心等式。第三種[V, D] eigs(A)注意這里函數(shù)名多了個(gè) s。eigs不是返回全部特征值而是按大小次序返回一部分通常用于大型稀疏矩陣。它底層用的是迭代法只算你指定的前 k 個(gè)。如果你對(duì)一個(gè) 10000×10000 的矩陣用eig它會(huì)試圖算全部特征值內(nèi)存和時(shí)間都可能爆掉用eigs指定算前 6 個(gè)就輕松很多。A [0 1 0 0; 1 0 0 0; 0 0 0 1; 0 0 1 0]; d eig(A); % 全部特征值列向量 [V, D] eig(A); % 特征向量矩陣 V 和特征值對(duì)角陣 D運(yùn)行后你會(huì)發(fā)現(xiàn)特征值是 1、-1、1、-1有重根。D 的主對(duì)角線上按順序排列V 的列向量對(duì)應(yīng)每個(gè)特征值。這里有個(gè)細(xì)節(jié)當(dāng)特征值有重根時(shí)eig返回的特征向量不一定是正交的甚至可能線性相關(guān)。如果你需要正交基得用實(shí)對(duì)稱矩陣那條路后面會(huì)講。2.3 結(jié)果怎么驗(yàn)證V^-1*A*V和符號(hào)化輸出算出[V, D]之后怎么確認(rèn)沒算錯(cuò)最直接的辦法是驗(yàn)算 V?1AV 是否等于 D。在 Matlab 里可以寫成V_inv inv(V); check V_inv * A * V; % 應(yīng)該接近對(duì)角陣 D如果 A 可對(duì)角化check的非對(duì)角元素應(yīng)該接近 0。但數(shù)值計(jì)算有浮點(diǎn)誤差你可能會(huì)看到 10?1? 量級(jí)的小數(shù)那是正常的。有時(shí)候你希望看到精確的分?jǐn)?shù)或根號(hào)形式而不是小數(shù)。這時(shí)可以把 V 轉(zhuǎn)成符號(hào)類型V_sym sym(V); V_sym^-1 * A * V_sym這樣輸出就是分?jǐn)?shù)和根號(hào)方便對(duì)照手算結(jié)果。但注意符號(hào)運(yùn)算只適合小矩陣驗(yàn)證大矩陣會(huì)非常慢。提示eig返回的特征向量默認(rèn)是歸一化的但歸一化方式不唯一。不同版本的 Matlab 可能在符號(hào)上差一個(gè)負(fù)號(hào)這不影響正確性。3. 對(duì)角化與若當(dāng)標(biāo)準(zhǔn)形能對(duì)角化最好不能也別慌3.1 可對(duì)角化的條件與實(shí)對(duì)稱矩陣的特殊待遇一個(gè) n 階方陣能不能對(duì)角化取決于它有沒有 n 個(gè)線性無關(guān)的特征向量。如果特征值都是單根那一定可對(duì)角化如果有重根就要看重根對(duì)應(yīng)的特征向量個(gè)數(shù)夠不夠。不夠的話就只能退而求其次求若當(dāng)標(biāo)準(zhǔn)形。實(shí)對(duì)稱矩陣是個(gè)特例它一定可以對(duì)角化而且特征向量可以選成兩兩正交的。Matlab 的eig對(duì)實(shí)對(duì)稱矩陣會(huì)自動(dòng)返回正交矩陣 P滿足 P?1AP D而且 P?1 就等于 P?。這意味著你不需要求逆轉(zhuǎn)置一下就行數(shù)值上更穩(wěn)定。A [2 2 -2; 2 5 -4; -2 -4 5]; % 實(shí)對(duì)稱矩陣 [P, D] eig(A); % P 是正交陣D 是對(duì)角陣這個(gè)性質(zhì)在二次型化簡里特別有用。二次型矩陣是實(shí)對(duì)稱的對(duì)角化之后就得到了標(biāo)準(zhǔn)形P 就是變換矩陣。你不需要手動(dòng)配方eig一步到位。3.2jordan函數(shù)與若當(dāng)塊的上三角定義當(dāng)矩陣不能對(duì)角化時(shí)jordan(A)返回它的若當(dāng)標(biāo)準(zhǔn)形 J。若當(dāng)標(biāo)準(zhǔn)形長什么樣對(duì)角線上是特征值對(duì)角線上方可能有 1形成一個(gè)個(gè)“塊”。Matlab 里若當(dāng)塊是按上三角形定義的也就是說 1 出現(xiàn)在對(duì)角線上方而不是下方。A [-1 -2 6; -1 0 3; -1 -1 4]; J jordan(A);運(yùn)行后 J 是一個(gè)上三角矩陣對(duì)角線上是特征值某些位置有 1。這些 1 的位置和個(gè)數(shù)反映了矩陣的“缺陷程度”——也就是特征向量不夠用的程度。但jordan有個(gè)血淚經(jīng)驗(yàn)它對(duì)數(shù)值誤差極其敏感。稍微擾動(dòng)一下矩陣元素若當(dāng)塊結(jié)構(gòu)就可能完全改變。所以實(shí)際工程里除非你確定矩陣是精確的符號(hào)矩陣或者小整數(shù)矩陣否則不要輕易依賴jordan的數(shù)值結(jié)果。它更適合教學(xué)和理論驗(yàn)證不適合處理實(shí)驗(yàn)測(cè)得的浮點(diǎn)數(shù)據(jù)。3.3 復(fù)對(duì)角矩陣轉(zhuǎn)實(shí)對(duì)角cdf2rdf在做什么實(shí)矩陣的特征值可能是復(fù)數(shù)而且成共軛對(duì)出現(xiàn)。eig返回的 D 可能是復(fù)對(duì)角陣對(duì)角線上有 abi 和 a?bi。如果你希望保持實(shí)數(shù)運(yùn)算可以用cdf2rdf把復(fù)對(duì)角陣轉(zhuǎn)成實(shí)塊對(duì)角陣。A [0 2 1; -2 0 3; -1 -3 0]; [v, d] eig(A); % d 是復(fù)對(duì)角陣 [V, D] cdf2rdf(v, d); % D 是實(shí)塊對(duì)角陣轉(zhuǎn)換之后D 的對(duì)角線上不再是復(fù)數(shù)而是 2×2 的實(shí)數(shù)塊。每個(gè) 2×2 塊對(duì)應(yīng)一對(duì)共軛特征值。這樣做的好處是后續(xù)計(jì)算全在實(shí)數(shù)域進(jìn)行避免復(fù)數(shù)運(yùn)算帶來的額外開銷和精度問題。注意cdf2rdf只對(duì)實(shí)矩陣有意義。如果原矩陣本身就是復(fù)矩陣這個(gè)轉(zhuǎn)換不適用。4. 常見誤區(qū)與易混概念別把這幾組函數(shù)搞混了4.1eig和eigs一個(gè)算全部一個(gè)算部分最容易混的就是這兩個(gè)。名字只差一個(gè)字母行為差很遠(yuǎn)。eig算全部特征值適合中小型稠密矩陣eigs算部分特征值適合大型稀疏矩陣。如果你對(duì)一個(gè) 5000×5000 的稀疏矩陣用eigMatlab 會(huì)試圖算 5000 個(gè)特征值內(nèi)存直接吃滿。用eigs(A, 6)只算前 6 個(gè)幾秒鐘就出結(jié)果。區(qū)分方法很簡單看矩陣大小和稀疏性。小矩陣、需要全部特征值用eig大矩陣、只關(guān)心最大或最小的幾個(gè)用eigs。4.2poly和polyfit名字像用途完全不同poly是求特征多項(xiàng)式系數(shù)的輸入是方陣輸出是系數(shù)向量。polyfit是多項(xiàng)式擬合輸入是數(shù)據(jù)點(diǎn)輸出是擬合系數(shù)。兩者都和“多項(xiàng)式”有關(guān)但一個(gè)服務(wù)于矩陣特征值一個(gè)服務(wù)于數(shù)據(jù)擬合。初學(xué)者看到poly容易聯(lián)想到擬合這是典型的望文生義。區(qū)分方法看輸入。輸入是方陣用poly輸入是 x 和 y 數(shù)據(jù)向量用polyfit。4.3jordan的數(shù)值陷阱理論成立不等于數(shù)值可靠理論上任何方陣都有若當(dāng)標(biāo)準(zhǔn)形但數(shù)值上jordan可能給你一個(gè)完全錯(cuò)誤的結(jié)果。原因在于若當(dāng)塊結(jié)構(gòu)對(duì)擾動(dòng)極其敏感。一個(gè)特征值接近重根的矩陣在浮點(diǎn)誤差下可能被判定為可對(duì)角化也可能被判定為有若當(dāng)塊結(jié)果不穩(wěn)定。區(qū)分方法如果矩陣元素是精確整數(shù)或符號(hào)jordan可信如果矩陣來自實(shí)驗(yàn)測(cè)量或浮點(diǎn)計(jì)算不要依賴jordan的結(jié)構(gòu)只看特征值就夠了。4.4balance不是“平衡方程”而是縮放矩陣balance(A)返回一個(gè)相似變換矩陣 T 和平衡矩陣 B使得 B T?1AT。它的作用是縮放矩陣的行和列讓各行各列的范數(shù)盡量接近從而提高特征值計(jì)算的精度。它不是用來解平衡方程的也不是用來做物理平衡的。區(qū)分方法看到balance出現(xiàn)在特征值計(jì)算上下文里就知道它在做矩陣縮放目的是數(shù)值穩(wěn)定。4.5trace和sum跡是特征值之和但sum不一定是trace(A)返回矩陣的跡等于主對(duì)角線元素之和也等于所有特征值之和。sum(A, dim)是元素求和dim1 按列求和dim2 按行求和。sum(sum(A,1),2)才是所有元素之和。初學(xué)者容易把trace和sum混為一談?dòng)绕涫窃隍?yàn)證特征值之和等于跡的時(shí)候。區(qū)分方法trace只針對(duì)方陣結(jié)果是標(biāo)量sum對(duì)任意矩陣都行結(jié)果可能是向量。5. 什么時(shí)候用、什么時(shí)候別用一個(gè)具體例子和邊界判斷假設(shè)你手頭有一個(gè) 3×3 的實(shí)對(duì)稱矩陣來自某個(gè)物理系統(tǒng)的剛度矩陣你想知道它的特征值和特征向量用來判斷系統(tǒng)的振動(dòng)模態(tài)。這時(shí)候eig是最合適的選擇因?yàn)閷?shí)對(duì)稱矩陣保證正交對(duì)角化數(shù)值穩(wěn)定結(jié)果可靠。A [2 2 -2; 2 5 -4; -2 -4 5]; [P, D] eig(A); % P 的列是振動(dòng)模態(tài)D 的對(duì)角線是頻率平方但如果你面對(duì)的是一個(gè) 20000×20000 的稀疏矩陣來自有限元分析你只關(guān)心前 10 個(gè)最小特征值。這時(shí)候絕對(duì)不能用eig必須用eigs(A, 10, sm)其中sm表示 smallest magnitude。否則內(nèi)存和時(shí)間都扛不住。再比如你有一個(gè)矩陣特征值有重根你想知道它能不能對(duì)角化。這時(shí)候不要直接上jordan因?yàn)閿?shù)值結(jié)果可能騙你。更穩(wěn)妥的做法是檢查特征向量的秩用[V, D] eig(A)然后看rank(V)是否等于矩陣階數(shù)。如果等于說明有足夠多的線性無關(guān)特征向量可對(duì)角化如果小于說明不能對(duì)角化若當(dāng)塊存在。還有一個(gè)場(chǎng)景你算完特征值發(fā)現(xiàn)結(jié)果里有很小的虛部比如 1.0000 0.0000i。這通常是浮點(diǎn)誤差不是真的復(fù)數(shù)特征值?？梢杂胷eal函數(shù)取實(shí)部或者用cdf2rdf轉(zhuǎn)成實(shí)塊對(duì)角陣。提示判斷特征值是否為實(shí)數(shù)不要看虛部是否嚴(yán)格為零而是看虛部是否小于某個(gè)容差比如abs(imag(d)) 1e-10。從那以后我每次用eig之前都會(huì)先問自己三個(gè)問題矩陣是稠密還是稀疏需要全部特征值還是前幾個(gè)矩陣有沒有特殊結(jié)構(gòu)對(duì)稱、三角、稀疏這三個(gè)問題決定了用eig還是eigs用不用balance以及結(jié)果怎么驗(yàn)證。希望幫到你。本文還有配套的精品資源點(diǎn)擊獲取