：特征值、若當標準形與符號運算避坑指南)
簡介這份資源面向?qū)W習MATLAB數(shù)值計算與線性代數(shù)的學生和工程人員聚焦矩陣處理中的特征值、特征向量與對角化等核心問題適合作為課堂實驗或自學練習的配套材料。壓縮包內(nèi)共1個doc文檔約98KB內(nèi)容以實驗講義形式組織涵蓋特征多項式求解、eig函數(shù)求特征值與特征向量、實對稱矩陣的正交對角化、balance函數(shù)提高計算精度以及jordan函數(shù)求若當標準形等知識點并配有可直接運行的MATLAB代碼示例與運算結(jié)果。文檔還整理了trace、sum、prod、cdf2rdf等常用矩陣函數(shù)的用法說明便于讀者對照實驗步驟復現(xiàn)結(jié)果、理解矩陣分析的基本流程。目前已有727人學習下載適合需要系統(tǒng)梳理MATLAB矩陣運算命令、完成相關(guān)實驗報告或查漏補缺的讀者參考使用。1. 從一次翻車的對角化說起這套 MATLAB 矩陣處理實驗到底值不值得跑帶過幾屆做數(shù)值計算方向的學生之后我發(fā)現(xiàn)一個規(guī)律凡是自稱學過線性代數(shù)的人第一次用 MATLAB 處理特征值問題十有八九會在對角化驗證那一步翻車。不是概念不懂而是從紙面推導到工程實現(xiàn)之間隔著一層東西——數(shù)值精度、符號運算、若當塊的排列方向每一個都能讓結(jié)果對不上。這份實驗六的資料恰好就是沖著這層隔閡來的。它圍繞特征值、特征向量、若當標準形這幾個矩陣分析里的硬骨頭把 MATLAB 里poly、eig、jordan、balance、cdf2rdf這些函數(shù)串成了一條完整的操作鏈。適合誰正在做線性代數(shù)實驗、需要把矩陣對角化落到代碼上的學生以及工作中偶爾要處理模態(tài)分析、系統(tǒng)穩(wěn)定性判據(jù)、二次型標準化的工程師。它不教你線性代數(shù)它教你怎么讓 MATLAB 把線性代數(shù)算對。2. 特征多項式與特征值poly、eig 和符號求解三條路怎么選2.1 特征多項式poly 的向量表示和 poly2str 的可讀化MATLAB 里求特征多項式最直接的是poly(A)它返回的是多項式系數(shù)的行向量按降冪排列。這個向量本身不好讀得配合poly2str轉(zhuǎn)成人類能看的字符串形式。% 定義方陣 A A [1 0; 2 3]; % poly 返回特征多項式的系數(shù)向量按降冪排列 p poly(A) % 把系數(shù)向量轉(zhuǎn)成可讀的多項式字符串 f poly2str(p, x)運行后p是[1 -4 3]對應x^2 - 4x 3。這里有個容易忽略的點poly(A)返回的系數(shù)向量長度是 n1首項系數(shù)恒為 1因為特征多項式是首一的。如果你拿這個向量去做roots得到的就是特征值。參數(shù)上沒什么可調(diào)的poly只接受方陣非方陣會直接報錯。另一種路子是從定義出發(fā)用符號變量硬算行列式A [1 0; 2 3]; E eye(2); % 2 階單位陣 syms x % 聲明符號變量 f det(x*E - A) % 直接按定義算特征多項式這條路的好處是結(jié)果以符號形式呈現(xiàn)(x-1)*(x-3)一眼就能看出根。代價是符號運算慢矩陣一大就卡。我一般只在矩陣階數(shù)小于 5 且需要精確表達式時才用符號法其余情況一律走poly加roots。2.2 eig 的四種調(diào)用形式與特征向量的歸一化問題eig是這套流程的核心。它至少有四種調(diào)用形式每種返回的東西不一樣混用是新手最常見的翻車點。A [0 1 0 0; 1 0 0 0; 0 0 0 1; 0 0 1 0]; % 形式一只返回特征值列向量 d eig(A) % 形式二返回特征向量矩陣 V 和特征值對角陣 D滿足 A*V V*D [V, D] eig(A) % 形式三eigs 按大小次序返回適合大型稀疏矩陣 [V, D] eigs(A) % 形式四廣義特征值問題 A*x lambda*B*x % [V, D] eig(A, B)形式二里有個關(guān)鍵細節(jié)D的主對角線元素是特征值但順序不保證按大小排。上面這個例子里D的對角線是-1, -1, 1, 1特征值 -1 是二重的對應V的第 1、2 列。這里V的列向量是歸一化的模長為 1但符號可能和你手算的不一樣——這是數(shù)值算法的正常行為不是錯誤。驗證對角化是否成立標準操作是% 用符號形式輸出 V避免浮點誤差干擾判斷 V_sym sym(V) % 驗證 V^-1 * A * V 是否等于 D V_sym^-1 * A * V_sym結(jié)果應該是一個對角陣對角線上就是特征值。如果這里算出來非對角元有微小數(shù)值別慌那是浮點誤差用sym轉(zhuǎn)成符號再算就干凈了。2.3 提高精度balance 什么時候該用balance這個函數(shù)很多人沒用過但它在特征值計算精度上確實能救場。原理是對矩陣做相似變換把行范數(shù)和列范數(shù)拉到接近的量級減少數(shù)值計算中的舍入誤差。% 返回相似變換矩陣 T 和平衡矩陣 B滿足 B T^-1 * A * T [T, B] balance(A); % 只要平衡矩陣 B balance(A);什么時候該用當你發(fā)現(xiàn)eig算出來的特征值對擾動特別敏感或者矩陣元素量級差異很大比如有的元素是 1e-8有的是 1e8先balance一下再求特征值精度會明顯改善。注意balance不改變特征值只改變矩陣的數(shù)值條件。我一般在對精度有要求且矩陣條件數(shù)偏大時才會加這一步日常小矩陣直接eig就夠了。3. 對角化實戰(zhàn)實對稱矩陣、復對角轉(zhuǎn)實對角與若當標準形3.1 實對稱矩陣的對角化eig 直接給出正交陣實對稱矩陣有個漂亮性質(zhì)一定可以對角化而且特征向量可以選成正交的。MATLAB 的eig對實對稱矩陣會自動返回正交矩陣P滿足P * A * P D。A [2 2 -2; 2 5 -4; -2 -4 5]; % 實對稱矩陣 [P, D] eig(A) % P 是正交陣D 是對角陣 % 驗證正交性P * P 應該是單位陣 P * P % 驗證對角化P * A * P 應該等于 D P * A * P運行后D的對角線是1, 1, 10P的列向量兩兩正交且模長為 1。這里P * P會得到一個接近單位陣的結(jié)果非對角元是 1e-16 量級的浮點噪聲。如果你在做二次型標準化P就是變換矩陣D就是標準形系數(shù)。注意實對稱矩陣不會出現(xiàn)復特征值所以不需要cdf2rdf處理。3.2 復對角矩陣轉(zhuǎn)實對角cdf2rdf 的 2x2 實數(shù)塊非對稱實矩陣的特征值可能是復數(shù)而且成共軛對出現(xiàn)。eig返回的D是對角陣但對角線上是復數(shù)。如果你需要實數(shù)形式的塊對角陣就得用cdf2rdf。A [0 2 1; -2 0 3; -1 -3 0]; [v, d] eig(A) % v 和 d 都是復數(shù)形式 [V, D] cdf2rdf(v, d) % 轉(zhuǎn)成實數(shù)塊對角形式轉(zhuǎn)換后D不再是對角陣而是在對角線上用 2x2 實數(shù)塊代替共軛復數(shù)對。上面這個例子里D會變成一個 3x3 的塊對角陣其中有一個 2x2 塊對應特征值0 ± 3.7417i。這個函數(shù)在控制系統(tǒng)里做模態(tài)分析時特別有用因為實矩陣形式更方便后續(xù)計算。參數(shù)上cdf2rdf要求輸入的v和d必須是eig輸出的復數(shù)形式不能先手動取實部。3.3 若當標準形jordan 函數(shù)與上三角塊的排列規(guī)則若當標準形是矩陣相似理論的終極形式MATLAB 里用jordan求。A [-1 -2 6; -1 0 3; -1 -1 4]; J jordan(A) % 返回 A 的若當標準形運行結(jié)果是J 1 1 0 0 1 0 0 0 1這里有個必須記住的坑MATLAB 的若當塊是按上三角形定義的也就是說非對角元 1 出現(xiàn)在對角元的上方。這和某些教材里下三角的定義相反。如果你拿教材上的結(jié)果直接對比發(fā)現(xiàn) 1 的位置不對先確認定義方向。另外jordan對數(shù)值矩陣的魯棒性一般矩陣接近虧損時結(jié)果可能不穩(wěn)定符號矩陣會好很多。3.4 跡與行列式的驗證trace、sum、prod 的配合特征值有兩個基本恒等式所有特征值之和等于矩陣的跡所有特征值之積等于矩陣的行列式。用 MATLAB 驗證這兩個性質(zhì)能快速判斷特征值算得對不對。A [1 2 2; 2 1 2; 2 2 1]; a trace(A) % 矩陣的跡 d eig(A) % 特征值 b sum(d, 1) % 特征值求和dim1 按列求和 e det(A) % 矩陣行列式 f prod(d, 1) % 特征值求積dim1 按列求積sum和prod的dim參數(shù)容易搞混dim1按列操作dim2按行操作。對于列向量dsum(d, 1)和sum(d, 2)結(jié)果一樣但習慣上求全部元素和用sum(sum(A, 1), 2)更穩(wěn)妥。驗證通過的標準是a和b相等、e和f相等浮點誤差在 1e-10 量級以內(nèi)都算正常。4. 避坑與排查特征值計算里那些讓人懷疑人生的報錯4.1 現(xiàn)象eig 返回的特征向量驗證不通過V^-1AV 不是對角陣原因通常有兩個。一是矩陣本身不可對角化虧損矩陣特征向量不夠 n 個線性無關(guān)的這時候V是奇異的V^-1根本不存在或者數(shù)值上極不穩(wěn)定。二是浮點誤差累積尤其是矩陣條件數(shù)很大時V^-1 * A * V的非對角元可能跑到 1e-8 量級看起來不像對角陣。解決辦法先用rank(V)檢查特征向量矩陣是否滿秩。如果不滿秩說明矩陣虧損得走若當標準形路線。如果是浮點誤差把V轉(zhuǎn)成符號形式sym(V)再算V^-1 * A * V結(jié)果就干凈了。我一般會在驗證前加一句format short控制顯示精度避免被一堆小數(shù)位干擾判斷。4.2 現(xiàn)象jordan 函數(shù)報錯或返回結(jié)果和教材對不上jordan對數(shù)值矩陣的容差很敏感。如果矩陣元素是浮點數(shù)且接近虧損MATLAB 可能判定為非虧損返回一個對角陣而不是帶若當塊的矩陣。另外前面提過若當塊的方向是上三角教材可能是下三角。解決辦法把矩陣轉(zhuǎn)成符號形式再調(diào)jordanjordan(sym(A))會走精確運算結(jié)果穩(wěn)定得多。方向問題只能靠對照定義解決記住 MATLAB 是上三角塊。如果矩陣規(guī)模大且必須用數(shù)值方法考慮先用balance改善條件數(shù)再算。4.3 現(xiàn)象cdf2rdf 報錯說輸入維度不匹配cdf2rdf要求輸入的v和d必須來自eig的復數(shù)輸出且d必須是對角陣形式。如果你手動構(gòu)造了一個實數(shù)矩陣傳進去或者d不是對角陣就會報錯。解決辦法確保調(diào)用鏈是[v, d] eig(A)之后直接[V, D] cdf2rdf(v, d)中間不要對v或d做任何修改。如果A是實對稱矩陣eig返回的本來就是實數(shù)不需要cdf2rdf強行調(diào)用反而會出錯。4.4 現(xiàn)象符號運算 solve 求特征值結(jié)果里有 root 表達式或復數(shù)用solve(f)求符號多項式的根時如果多項式次數(shù)大于 4MATLAB 可能返回root形式的符號解或者帶復數(shù)。這不是錯誤是代數(shù)基本定理的必然結(jié)果——五次以上沒有通用根式解。解決辦法對高次多項式直接用roots求數(shù)值解或者用vpa把符號解轉(zhuǎn)成高精度數(shù)值。如果只關(guān)心實特征值用real取實部后再判斷虛部是否可忽略。我一般對超過 4 階的矩陣直接走eig不碰符號求解省得跟root表達式較勁。4.5 現(xiàn)象balance 之后特征值變了balance做的是相似變換理論上不改變特征值。如果你發(fā)現(xiàn)變了大概率是變換后的矩陣條件數(shù)改善數(shù)值計算更準了所以結(jié)果和沒平衡之前有微小差異。這不是 bug是精度提升的表現(xiàn)。解決辦法把balance前后的特征值都算出來對比差異應該在合理范圍內(nèi)。如果差異巨大檢查T矩陣是否可逆balance返回的T在極端情況下可能接近奇異。我一般只在確認矩陣條件數(shù)很差時才用balance日常計算不加這一步。5. 進階技巧用符號工具箱做精確對角化與二次型標準化前面幾章走的是數(shù)值路線快但有浮點誤差。如果你需要精確結(jié)果——比如做二次型標準化、推導變換矩陣的解析形式——符號工具箱是更好的選擇。下面這套流程是我在需要干凈結(jié)果時的固定操作。先看一個完整的符號對角化例子A [1 2 2; 2 1 2; 2 2 1]; E eye(3); syms x f det(x*E - A) % 特征多項式符號形式 eigenvalues solve(f) % 精確特征值f會返回x^3 - 3*x^2 - 9*x - 5solve給出5, -1, -1。接下來對每個特征值求特征向量% 特征值 5 對應的特征向量 y 5; B y*E - A; b1 sym(null(B)) % 基礎解系符號形式 % 特征值 -1 對應的特征向量 y -1; B y*E - A; b2 sym(null(B)) % 返回兩個基向量null(B)返回的是齊次方程組B*Z 0的基礎解系。對于二重特征值 -1b2是一個 3x2 的矩陣兩列分別是兩個線性無關(guān)的特征向量。把b1和b2拼起來得到變換矩陣TT [b1, b2] % 所有特征向量拼成的矩陣 D T^-1 * A * T % 驗證對角化結(jié)果應該是精確對角陣D會返回一個精確的對角陣對角線上是5, -1, -1非對角元嚴格為零。這就是符號運算的價值——沒有浮點噪聲結(jié)果可以直接寫進推導過程。再看二次型標準化的場景。給定二次型矩陣A要把它化成標準形步驟是步驟操作對應代碼1確認矩陣對稱isequal(A, A)2求特征值和正交特征向量[P, D] eig(A)3構(gòu)造正交變換x P * y4代入二次型得到標準形y * D * y5驗證變換矩陣正交性P * P接近單位陣這套流程的關(guān)鍵在第 2 步實對稱矩陣的eig返回的P一定是正交陣所以P^-1就是P不用真的求逆。如果你拿到的P不正交檢查矩陣是否真的對稱——數(shù)值上A - A的范數(shù)應該在 1e-12 以下才算對稱。最后一個技巧是關(guān)于精度控制的。當你需要高精度數(shù)值結(jié)果時vpa配合eig的符號輸出能派上用場A [1 2 2; 2 1 2; 2 2 1]; [V, D] eig(sym(A)); % 符號形式求特征值 V_num double(V); % 轉(zhuǎn)回雙精度 D_num double(D); % 需要更高精度時用 vpa V_high vpa(V, 32); % 32 位有效數(shù)字vpa的第二個參數(shù)控制有效數(shù)字位數(shù)默認是 32。位數(shù)越高計算越慢一般 32 位足夠工程使用。我通常在驗證理論推導時才用vpa日常計算double就夠了。從那以后我每次做對角化驗證都強制走一遍符號轉(zhuǎn)數(shù)值再對比的流程——先用sym算精確結(jié)果再用double轉(zhuǎn)回來對比數(shù)值解差異超過 1e-10 就回頭查矩陣條件數(shù)。這個習慣幫我逮住過好幾次隱藏的虧損矩陣問題。希望幫到你。本文還有配套的精品資源點擊獲取 MATLAB矩陣對角化實戰(zhàn)：特征值、若當標準形與符號運算避坑指南
簡介這份資源面向?qū)W習MATLAB數(shù)值計算與線性代數(shù)的學生和工程人員聚焦矩陣處理中的特征值、特征向量與對角化等核心問題適合作為課堂實驗或自學練習的配套材料。壓縮包內(nèi)共1個doc文檔約98KB內(nèi)容以實驗講義形式組織涵蓋特征多項式求解、eig函數(shù)求特征值與特征向量、實對稱矩陣的正交對角化、balance函數(shù)提高計算精度以及jordan函數(shù)求若當標準形等知識點并配有可直接運行的MATLAB代碼示例與運算結(jié)果。文檔還整理了trace、sum、prod、cdf2rdf等常用矩陣函數(shù)的用法說明便于讀者對照實驗步驟復現(xiàn)結(jié)果、理解矩陣分析的基本流程。目前已有727人學習下載適合需要系統(tǒng)梳理MATLAB矩陣運算命令、完成相關(guān)實驗報告或查漏補缺的讀者參考使用。1. 從一次翻車的對角化說起這套 MATLAB 矩陣處理實驗到底值不值得跑帶過幾屆做數(shù)值計算方向的學生之后我發(fā)現(xiàn)一個規(guī)律凡是自稱學過線性代數(shù)的人第一次用 MATLAB 處理特征值問題十有八九會在對角化驗證那一步翻車。不是概念不懂而是從紙面推導到工程實現(xiàn)之間隔著一層東西——數(shù)值精度、符號運算、若當塊的排列方向每一個都能讓結(jié)果對不上。這份實驗六的資料恰好就是沖著這層隔閡來的。它圍繞特征值、特征向量、若當標準形這幾個矩陣分析里的硬骨頭把 MATLAB 里poly、eig、jordan、balance、cdf2rdf這些函數(shù)串成了一條完整的操作鏈。適合誰正在做線性代數(shù)實驗、需要把矩陣對角化落到代碼上的學生以及工作中偶爾要處理模態(tài)分析、系統(tǒng)穩(wěn)定性判據(jù)、二次型標準化的工程師。它不教你線性代數(shù)它教你怎么讓 MATLAB 把線性代數(shù)算對。2. 特征多項式與特征值poly、eig 和符號求解三條路怎么選2.1 特征多項式poly 的向量表示和 poly2str 的可讀化MATLAB 里求特征多項式最直接的是poly(A)它返回的是多項式系數(shù)的行向量按降冪排列。這個向量本身不好讀得配合poly2str轉(zhuǎn)成人類能看的字符串形式。% 定義方陣 A A [1 0; 2 3]; % poly 返回特征多項式的系數(shù)向量按降冪排列 p poly(A) % 把系數(shù)向量轉(zhuǎn)成可讀的多項式字符串 f poly2str(p, x)運行后p是[1 -4 3]對應x^2 - 4x 3。這里有個容易忽略的點poly(A)返回的系數(shù)向量長度是 n1首項系數(shù)恒為 1因為特征多項式是首一的。如果你拿這個向量去做roots得到的就是特征值。參數(shù)上沒什么可調(diào)的poly只接受方陣非方陣會直接報錯。另一種路子是從定義出發(fā)用符號變量硬算行列式A [1 0; 2 3]; E eye(2); % 2 階單位陣 syms x % 聲明符號變量 f det(x*E - A) % 直接按定義算特征多項式這條路的好處是結(jié)果以符號形式呈現(xiàn)(x-1)*(x-3)一眼就能看出根。代價是符號運算慢矩陣一大就卡。我一般只在矩陣階數(shù)小于 5 且需要精確表達式時才用符號法其余情況一律走poly加roots。2.2 eig 的四種調(diào)用形式與特征向量的歸一化問題eig是這套流程的核心。它至少有四種調(diào)用形式每種返回的東西不一樣混用是新手最常見的翻車點。A [0 1 0 0; 1 0 0 0; 0 0 0 1; 0 0 1 0]; % 形式一只返回特征值列向量 d eig(A) % 形式二返回特征向量矩陣 V 和特征值對角陣 D滿足 A*V V*D [V, D] eig(A) % 形式三eigs 按大小次序返回適合大型稀疏矩陣 [V, D] eigs(A) % 形式四廣義特征值問題 A*x lambda*B*x % [V, D] eig(A, B)形式二里有個關(guān)鍵細節(jié)D的主對角線元素是特征值但順序不保證按大小排。上面這個例子里D的對角線是-1, -1, 1, 1特征值 -1 是二重的對應V的第 1、2 列。這里V的列向量是歸一化的模長為 1但符號可能和你手算的不一樣——這是數(shù)值算法的正常行為不是錯誤。驗證對角化是否成立標準操作是% 用符號形式輸出 V避免浮點誤差干擾判斷 V_sym sym(V) % 驗證 V^-1 * A * V 是否等于 D V_sym^-1 * A * V_sym結(jié)果應該是一個對角陣對角線上就是特征值。如果這里算出來非對角元有微小數(shù)值別慌那是浮點誤差用sym轉(zhuǎn)成符號再算就干凈了。2.3 提高精度balance 什么時候該用balance這個函數(shù)很多人沒用過但它在特征值計算精度上確實能救場。原理是對矩陣做相似變換把行范數(shù)和列范數(shù)拉到接近的量級減少數(shù)值計算中的舍入誤差。% 返回相似變換矩陣 T 和平衡矩陣 B滿足 B T^-1 * A * T [T, B] balance(A); % 只要平衡矩陣 B balance(A);什么時候該用當你發(fā)現(xiàn)eig算出來的特征值對擾動特別敏感或者矩陣元素量級差異很大比如有的元素是 1e-8有的是 1e8先balance一下再求特征值精度會明顯改善。注意balance不改變特征值只改變矩陣的數(shù)值條件。我一般在對精度有要求且矩陣條件數(shù)偏大時才會加這一步日常小矩陣直接eig就夠了。3. 對角化實戰(zhàn)實對稱矩陣、復對角轉(zhuǎn)實對角與若當標準形3.1 實對稱矩陣的對角化eig 直接給出正交陣實對稱矩陣有個漂亮性質(zhì)一定可以對角化而且特征向量可以選成正交的。MATLAB 的eig對實對稱矩陣會自動返回正交矩陣P滿足P * A * P D。A [2 2 -2; 2 5 -4; -2 -4 5]; % 實對稱矩陣 [P, D] eig(A) % P 是正交陣D 是對角陣 % 驗證正交性P * P 應該是單位陣 P * P % 驗證對角化P * A * P 應該等于 D P * A * P運行后D的對角線是1, 1, 10P的列向量兩兩正交且模長為 1。這里P * P會得到一個接近單位陣的結(jié)果非對角元是 1e-16 量級的浮點噪聲。如果你在做二次型標準化P就是變換矩陣D就是標準形系數(shù)。注意實對稱矩陣不會出現(xiàn)復特征值所以不需要cdf2rdf處理。3.2 復對角矩陣轉(zhuǎn)實對角cdf2rdf 的 2x2 實數(shù)塊非對稱實矩陣的特征值可能是復數(shù)而且成共軛對出現(xiàn)。eig返回的D是對角陣但對角線上是復數(shù)。如果你需要實數(shù)形式的塊對角陣就得用cdf2rdf。A [0 2 1; -2 0 3; -1 -3 0]; [v, d] eig(A) % v 和 d 都是復數(shù)形式 [V, D] cdf2rdf(v, d) % 轉(zhuǎn)成實數(shù)塊對角形式轉(zhuǎn)換后D不再是對角陣而是在對角線上用 2x2 實數(shù)塊代替共軛復數(shù)對。上面這個例子里D會變成一個 3x3 的塊對角陣其中有一個 2x2 塊對應特征值0 ± 3.7417i。這個函數(shù)在控制系統(tǒng)里做模態(tài)分析時特別有用因為實矩陣形式更方便后續(xù)計算。參數(shù)上cdf2rdf要求輸入的v和d必須是eig輸出的復數(shù)形式不能先手動取實部。3.3 若當標準形jordan 函數(shù)與上三角塊的排列規(guī)則若當標準形是矩陣相似理論的終極形式MATLAB 里用jordan求。A [-1 -2 6; -1 0 3; -1 -1 4]; J jordan(A) % 返回 A 的若當標準形運行結(jié)果是J 1 1 0 0 1 0 0 0 1這里有個必須記住的坑MATLAB 的若當塊是按上三角形定義的也就是說非對角元 1 出現(xiàn)在對角元的上方。這和某些教材里下三角的定義相反。如果你拿教材上的結(jié)果直接對比發(fā)現(xiàn) 1 的位置不對先確認定義方向。另外jordan對數(shù)值矩陣的魯棒性一般矩陣接近虧損時結(jié)果可能不穩(wěn)定符號矩陣會好很多。3.4 跡與行列式的驗證trace、sum、prod 的配合特征值有兩個基本恒等式所有特征值之和等于矩陣的跡所有特征值之積等于矩陣的行列式。用 MATLAB 驗證這兩個性質(zhì)能快速判斷特征值算得對不對。A [1 2 2; 2 1 2; 2 2 1]; a trace(A) % 矩陣的跡 d eig(A) % 特征值 b sum(d, 1) % 特征值求和dim1 按列求和 e det(A) % 矩陣行列式 f prod(d, 1) % 特征值求積dim1 按列求積sum和prod的dim參數(shù)容易搞混dim1按列操作dim2按行操作。對于列向量dsum(d, 1)和sum(d, 2)結(jié)果一樣但習慣上求全部元素和用sum(sum(A, 1), 2)更穩(wěn)妥。驗證通過的標準是a和b相等、e和f相等浮點誤差在 1e-10 量級以內(nèi)都算正常。4. 避坑與排查特征值計算里那些讓人懷疑人生的報錯4.1 現(xiàn)象eig 返回的特征向量驗證不通過V^-1AV 不是對角陣原因通常有兩個。一是矩陣本身不可對角化虧損矩陣特征向量不夠 n 個線性無關(guān)的這時候V是奇異的V^-1根本不存在或者數(shù)值上極不穩(wěn)定。二是浮點誤差累積尤其是矩陣條件數(shù)很大時V^-1 * A * V的非對角元可能跑到 1e-8 量級看起來不像對角陣。解決辦法先用rank(V)檢查特征向量矩陣是否滿秩。如果不滿秩說明矩陣虧損得走若當標準形路線。如果是浮點誤差把V轉(zhuǎn)成符號形式sym(V)再算V^-1 * A * V結(jié)果就干凈了。我一般會在驗證前加一句format short控制顯示精度避免被一堆小數(shù)位干擾判斷。4.2 現(xiàn)象jordan 函數(shù)報錯或返回結(jié)果和教材對不上jordan對數(shù)值矩陣的容差很敏感。如果矩陣元素是浮點數(shù)且接近虧損MATLAB 可能判定為非虧損返回一個對角陣而不是帶若當塊的矩陣。另外前面提過若當塊的方向是上三角教材可能是下三角。解決辦法把矩陣轉(zhuǎn)成符號形式再調(diào)jordanjordan(sym(A))會走精確運算結(jié)果穩(wěn)定得多。方向問題只能靠對照定義解決記住 MATLAB 是上三角塊。如果矩陣規(guī)模大且必須用數(shù)值方法考慮先用balance改善條件數(shù)再算。4.3 現(xiàn)象cdf2rdf 報錯說輸入維度不匹配cdf2rdf要求輸入的v和d必須來自eig的復數(shù)輸出且d必須是對角陣形式。如果你手動構(gòu)造了一個實數(shù)矩陣傳進去或者d不是對角陣就會報錯。解決辦法確保調(diào)用鏈是[v, d] eig(A)之后直接[V, D] cdf2rdf(v, d)中間不要對v或d做任何修改。如果A是實對稱矩陣eig返回的本來就是實數(shù)不需要cdf2rdf強行調(diào)用反而會出錯。4.4 現(xiàn)象符號運算 solve 求特征值結(jié)果里有 root 表達式或復數(shù)用solve(f)求符號多項式的根時如果多項式次數(shù)大于 4MATLAB 可能返回root形式的符號解或者帶復數(shù)。這不是錯誤是代數(shù)基本定理的必然結(jié)果——五次以上沒有通用根式解。解決辦法對高次多項式直接用roots求數(shù)值解或者用vpa把符號解轉(zhuǎn)成高精度數(shù)值。如果只關(guān)心實特征值用real取實部后再判斷虛部是否可忽略。我一般對超過 4 階的矩陣直接走eig不碰符號求解省得跟root表達式較勁。4.5 現(xiàn)象balance 之后特征值變了balance做的是相似變換理論上不改變特征值。如果你發(fā)現(xiàn)變了大概率是變換后的矩陣條件數(shù)改善數(shù)值計算更準了所以結(jié)果和沒平衡之前有微小差異。這不是 bug是精度提升的表現(xiàn)。解決辦法把balance前后的特征值都算出來對比差異應該在合理范圍內(nèi)。如果差異巨大檢查T矩陣是否可逆balance返回的T在極端情況下可能接近奇異。我一般只在確認矩陣條件數(shù)很差時才用balance日常計算不加這一步。5. 進階技巧用符號工具箱做精確對角化與二次型標準化前面幾章走的是數(shù)值路線快但有浮點誤差。如果你需要精確結(jié)果——比如做二次型標準化、推導變換矩陣的解析形式——符號工具箱是更好的選擇。下面這套流程是我在需要干凈結(jié)果時的固定操作。先看一個完整的符號對角化例子A [1 2 2; 2 1 2; 2 2 1]; E eye(3); syms x f det(x*E - A) % 特征多項式符號形式 eigenvalues solve(f) % 精確特征值f會返回x^3 - 3*x^2 - 9*x - 5solve給出5, -1, -1。接下來對每個特征值求特征向量% 特征值 5 對應的特征向量 y 5; B y*E - A; b1 sym(null(B)) % 基礎解系符號形式 % 特征值 -1 對應的特征向量 y -1; B y*E - A; b2 sym(null(B)) % 返回兩個基向量null(B)返回的是齊次方程組B*Z 0的基礎解系。對于二重特征值 -1b2是一個 3x2 的矩陣兩列分別是兩個線性無關(guān)的特征向量。把b1和b2拼起來得到變換矩陣TT [b1, b2] % 所有特征向量拼成的矩陣 D T^-1 * A * T % 驗證對角化結(jié)果應該是精確對角陣D會返回一個精確的對角陣對角線上是5, -1, -1非對角元嚴格為零。這就是符號運算的價值——沒有浮點噪聲結(jié)果可以直接寫進推導過程。再看二次型標準化的場景。給定二次型矩陣A要把它化成標準形步驟是步驟操作對應代碼1確認矩陣對稱isequal(A, A)2求特征值和正交特征向量[P, D] eig(A)3構(gòu)造正交變換x P * y4代入二次型得到標準形y * D * y5驗證變換矩陣正交性P * P接近單位陣這套流程的關(guān)鍵在第 2 步實對稱矩陣的eig返回的P一定是正交陣所以P^-1就是P不用真的求逆。如果你拿到的P不正交檢查矩陣是否真的對稱——數(shù)值上A - A的范數(shù)應該在 1e-12 以下才算對稱。最后一個技巧是關(guān)于精度控制的。當你需要高精度數(shù)值結(jié)果時vpa配合eig的符號輸出能派上用場A [1 2 2; 2 1 2; 2 2 1]; [V, D] eig(sym(A)); % 符號形式求特征值 V_num double(V); % 轉(zhuǎn)回雙精度 D_num double(D); % 需要更高精度時用 vpa V_high vpa(V, 32); % 32 位有效數(shù)字vpa的第二個參數(shù)控制有效數(shù)字位數(shù)默認是 32。位數(shù)越高計算越慢一般 32 位足夠工程使用。我通常在驗證理論推導時才用vpa日常計算double就夠了。從那以后我每次做對角化驗證都強制走一遍符號轉(zhuǎn)數(shù)值再對比的流程——先用sym算精確結(jié)果再用double轉(zhuǎn)回來對比數(shù)值解差異超過 1e-10 就回頭查矩陣條件數(shù)。這個習慣幫我逮住過好幾次隱藏的虧損矩陣問題。希望幫到你。本文還有配套的精品資源點擊獲取