1. 程式人生 > >矩陣的五種分解的matlab實現

矩陣的五種分解的matlab實現

由於這學期修了矩陣分析這門課,課程要求用matlab實現矩陣的5種分解,僅僅是實現了分解,上傳到部落格存檔,萬一哪天某位同學就需要了呢。。

 

1.矩陣的滿秩分解

  • 程式碼實現
 1 %矩陣的滿秩分解
 2 clear
 3 %設輸入矩陣為M(P152 例4.1.1 4 A = [1,4,-1,5,6;
 5     2,0,0,0,-14;
 6     -1,2,-4,0,1;
 7     2,6,-5,5,-7]
 8 A1 = rref(A);    %將矩陣A化成行最簡形式儲存在A1中
 9 [m,n]=size(A);    %獲取矩陣A的大小:m行n列
10 B0= [];%生成一個空向量 11 C0= [];%生成一個空向量 12 for i=1:m %依次掃描矩陣m行 13 flag=1; 14 for j=1:n %依次掃描矩陣n列 15 if A1(i,j)==1 %若A1(i, j)等於1 16 for k=1:i-1 %固定j列,掃描此列的第1行到i-1行元素 17 if A1(k,j)~=0 %判斷是否全為0 18 flag=0; %若不全為0,則將flag置為0(說明此列不是單位矩陣的列)
19 break; 20 end 21 end 22 for k=i+1:m %固定j列,掃描此列的第i+1行到m行(即最後一行)元素 23 if A1(k,j)~=0 %判斷是否全為0 24 flag=0; %若不全為0,則將flag置為0(說明此列不是單位矩陣的列) 25 break; 26 end 27 end
28 if flag==1 %若flag為1(不為0),則說明此列是【矩陣的行最簡形式矩陣】的單位矩陣的列 29 B0=[B0,A(:,j)]; %將矩陣A的j列加到B0列向量之後 30 C0=[C0;A1(i,:)]; %將矩陣A1的i行加到C0行向量之後, 31 end 32 end 33 end 34 end 35 [m1,n1]=size(B0); %獲取矩陣B0的大小:m1行n1列 36 [m2,n2]=size(C0); %獲取矩陣C0的大小:m2行n2列 37 B=B0(:,1:n1) %將矩陣B0的第1列到最後一列賦值給矩陣B 38 C=C0(1:m2,:) %將矩陣C0的第1行到最後一行賦值給矩陣C 39 %驗證:BC=A 40 A_1= B*C

 

2.矩陣的正交三角分解

  • 程式碼實現

直接呼叫matlab自帶qr()函式即可

1 %矩陣的正交三角分解
2 clear;
3 A = [-3,1,-2;1,1,1;1,-1,0;1,-1,1]
4 [Q, R] = qr(A) %正交三角分解,Q為酉矩陣,R為正交下三角矩陣
5 %驗證:QR是否為A,以及Q是否為酉矩陣
6 A_1 = Q * R
7 Q_1 = Q * conj(Q.')

 

3.矩陣的奇異值分解

  • 程式碼實現
1 %矩陣的奇異值分解
2 clear,clc
3 A = [1,1;0,0;1,1];
4 [U,S,V] = svd(A) %返回一個與A同大小的對角矩陣S,兩個酉矩陣U和V,且滿足A= U*S*V~H。
5                  %若A為m×n陣,則U為m×m陣,V為n×n陣。奇異值在S的對角線上,非負且按降序排列。
6 
7 %驗證A=USV~H
8 A = [1,1;0,0;1,1]
9 A_1 = U*S*conj(V.') 

 

4.矩陣的極分解

  • 程式碼實現
1 %矩陣的極分解
2 clear,clc;
3 A = [2,1,2;0,1,3;1,0,0];
4 H1 = sqrtm(A*A') %返回矩陣的主要平方根
5 U1 = inv(H1)*A %求逆
6 A_1 = H1*U1 
7 H2 = sqrtm(A)
8 U2 = A*inv(H2)
9 A_2 = U2*H2

 

5.矩陣的譜分解

以正規矩陣為例:

 

  •  程式碼實現
 1 %矩陣的譜分解
 2 clear,clc
 3 A = [4,6,0;-3,-5,0;-3,-6,1]; %單純矩陣
 4 %A = [-2i,4,-2;-4,-3i,-3i;2,-2i,-5i]; %正規矩陣
 5 [V,D] = eig(A) %求特徵值與特徵向量
 6 
 7 %正交歸一化
 8 V_C = orth(V) ;%特徵向量正交化
 9 V_C_Z = V_C./repmat(sqrt(sum(V_C.^2,1)),size(V_C,1),1); %特徵向量列歸一化
10 
11 A_H = A * conj(A');%求A的共軛轉置
12 if A == A_H | A == -(A_H)  %判斷是否是正規矩陣
13     [m,n] = size(V_C_Z);
14     G2 = zeros(m,n); 
15     for i=1:n
16         G1 = V_C_Z(:,i) * conj(V_C_Z');
17         G2 = G2 + D(i, i) * G1;
18     end
19     G_Z = G2
20 else  %否則是單純矩陣
21     P_1 = (inv(V))';
22     [m,n] = size(P_1);
23     G3 = zeros(m,n);
24     for i=1:n 
25         G4 = V(:,i) * (P_1(:,i))';
26         G3 = G3 + D(i,i)* G4;
27     end 
28     G_R = G3
29 end

 

 PS:滿秩分解的參考地址記不住了,這裡就不備註了,僅僅出於學習的目的,不喜勿噴。