張智星 (Roger Jang) 台大資工系 多媒體檢索實驗室

Presentation on theme: "張智星 (Roger Jang) 台大資工系 多媒體檢索實驗室"— Presentation transcript:

1 張智星 (Roger Jang) jang@mirlab.org http://mirlab.org/jang 台大資工系 多媒體檢索實驗室
MATLAB 程式設計進階篇 稀疏矩陣 張智星 (Roger Jang) 台大資工系 多媒體檢索實驗室

2 5-1 稀疏矩陣的建立 根據內部資料結構的不同，MATLAB 矩陣分為兩種 : 完全矩陣（Full Matrix）
5-1 稀疏矩陣的建立 根據內部資料結構的不同，MATLAB 矩陣分為兩種 : 完全矩陣（Full Matrix） 每個元素都有自己的儲存空間 每個元素的佔用空間大小依其資料型態而有所不同 double8 bytes, int324 bytes, int162 bytes, int81 bytes 一個 m×n 的完全矩陣，所佔用的記憶體空間是 k×m×n bytes 稀疏矩陣（Sparse Matrix） 大部分的元素都是 0 只須儲存「非零元素的位置」及其「元素值」 優點 節省記憶體儲存空間稀疏矩陣的最重要任務！ 節省某一些不必要的運算

3 稀疏矩陣範例-1 (I) sparse 指令可將一個完全矩陣轉換成稀疏矩陣 範例5-1：sparse01.m
(1,1) (3,2) (2,4) 可看出，S 是一個稀疏矩陣，MATLAB 只儲存其各個非零元素的位置（即其二維下標 (1,1)、(3,2)、(2,4)）和元素值（2、4、1） clear all; % 清除所有的變數 A = [ ; ; ]; % 完全矩陣 S = sparse(A) % 將完全矩陣 A 轉換成稀疏矩陣 S

4 稀疏矩陣範例-1 (II) 比較矩陣 A 和 S 所佔用的記憶體大小 :
>> whos Name Size Bytes Class Attributes A x double S x double sparse 看出稀疏矩陣 S 佔用記憶體的位元組數目（88 Bytes）比 矩陣A（佔用 96 Bytes）還要小。 當矩陣變大、密度變低時，這個差異會越來越大。

5 稀疏矩陣範例-2(I) S = sparse(i, j, s, m, n); 使用 sparse 指令來直接產生稀疏矩陣 i 是列索引
m 是 s 的列維度 n 是 s 的行維度 i、j、s 都是長度相同的向量 s(k) 的二維下標即是 i(k)及 j(k)

6 稀疏矩陣範例-2(II) >> S = sparse([1 3 2], [1 2 4], [2 4 1], 3, 4) S = (1,1) (3,2) (2,4) 也可以在 sparse 指令加上第六種輸入變數，代表最多可以容納的非零元素個素，使得您可以後續再加入非零元素，而不必改變整個稀疏矩陣的結構

7 稀疏矩陣範例-3(I) 指令 spdiags，可由對角線元素來建構一個稀疏矩陣 S = spdiags(D, p, m, n)

8 稀疏矩陣範例-3(II) 範例5-2：sparse02.m S = D = [3 2 9; 2 4 9; 1 1 4];
(1,1) (2,1) (2,2) (3,2) (1,3) (3,3) (4,3) D = [3 2 9; 2 4 9; 1 1 4]; d = [2 0 -1]; S = spdiags(D, d, 4, 3)

9 稀疏矩陣範例-3(III) 以完全矩陣來顯示，可得 : >> A = full(S)
可看出矩陣 A 的三個非零對角線向量分別是 D 的三個行向量。（在D的第一個行向量中，只用到最後一個元素1。） 提示： D = d =

10 稀疏矩陣範例-4(I) 一般的 load 及 save 指令，也可以處理稀疏矩陣，並儲存於二進制的 MAT 檔案
spconvert指令，可將一個 m×3 的矩陣轉換成稀疏矩陣 第一直行代表列索引 第二直行代表行索引 第三直行則是非零的元素值 檔案 spmat.dat 的內容可顯示如下： >> type spmat.dat 1 1 2 3 2 4 2 4 1 8 3 5

11 稀疏矩陣範例-4(II) 建構此稀疏矩陣 範例5-3：spconvert01.m S = load spmat.dat (1,1) 2
(1,1) (3,2) (8,3) (2,4) load spmat.dat S = spconvert(spmat)

12 5-2 稀疏矩陣的儲存空間 一個只包含實數的稀疏矩陣，假設其維度為 m×n，含有 nnz 個非零元素，MATLAB 動用了三個內部陣列來儲存此稀疏矩陣的相關資訊： 第一個陣列：以 double 方式儲存了所有的非零元素，其長度為 nnz，使用的空間為大小 8*nnz 位元組（Bytes）。 第二個陣列：以整數方式儲存了每個元素的列索引，其長度為 nnz，使用的空間大小為 4*nnz 位元組（Bytes）。 第三個陣列：以整數方式儲存了每個直行的起始指標，其長度為 n，使用空間大小為 4*n 位元組（Bytes）。 注意：上述儲存方式僅是用於MATLAB第六版！

13 稀疏矩陣的儲存空間 MATLAB第七版後，整個稀疏矩陣佔用的空間大小為 16*nnz + 8*n + 8，可以驗證如下。
範例5-4：memorySize01.m totalSize = 34032 Name Size Bytes Class Attributes west x double sparse clear all load west0479.mat [m,n]=size(west0479); totalSize=16*nnz(west0479)+8*n+8 whos west0479

14 稀疏矩陣的儲存空間 Regarding the formula “16*nnz + 8*n + 8”
Can be verified via least-square estimate in memorySize02.m Open question: Why is it so?

15 稀疏矩陣的相關指令 其他相關指令： nnz：可傳回稀疏矩陣的非零元素個數 nonzeros：傳回一個包含所有非零元素的行向量
nzmax：傳回最大的非零元素個數

16 提示 在一個稀疏矩陣中，將一個非零元素設定成零時，MATLAB 並不會自動釋放記憶體空間，換包話說，nnz 會隨著非零元素的減少而減少，但 nzmax 並不會隨著改變 但是多加一個非零元素時，若 nnz 已大於 nzmax 時，nzmax 會隨之增大（即 MATLAB 會自動配置記憶體）以儲存新加的元素

17 5-3 稀疏矩陣的觀看與圖示 spy 指令可用於觀看稀疏矩陣的非零元素分佈情況 範例5-5：spy01.m
5-3 稀疏矩陣的觀看與圖示 spy 指令可用於觀看稀疏矩陣的非零元素分佈情況 範例5-5：spy01.m load west0479.mat % 載入二進位制檔案 west0479.mat spy(west0479) % 觀看稀疏矩陣的非零元素分佈情況

18 稀疏矩陣的觀看與圖示 矩陣 west0479 的維度是 479*479，但是只包含 1887 個非零元素，因此此矩陣的密度只有 1887/(479*479) =

19 稀疏矩陣的觀看與圖示 稀疏矩陣特別適用於表示一個「無向圖」（Undirected Graph）的「鄰近矩陣」（Adjacency Matrix） 若某圖的第 i 和第 j 個節點有直線連接，則其相對應的鄰近矩陣在第 i 列、第 j 行的元素值為 1，其他元素值則為零 為節省空間，僅用上三角或是下三角矩陣來表示

20 稀疏矩陣的觀看與圖示 對應的鄰近矩陣可表示成：
>> A = spconvert([1 2 1; 2 3 1; 2 4 1; 3 4 1; 3 5 1; 5 6 1; 4 6 1]); >> full(A) ans =

21 稀疏矩陣的觀看與圖示 假設這 6 個節點的座標如下： 可用 gplot 指令來畫出上述的無向圖： 範例5-6：gplot01.m
>> xy = [0 1; 1 2; 1 0; 2 0; 2 2; 3 1]; % 每個列向量是一組座標 可用 gplot 指令來畫出上述的無向圖： 範例5-6：gplot01.m 其中 '-o' 代表以實線（'-'）及圓圈（'o'）來作圖 A = spconvert([1 2 1; 2 3 1; 2 4 1; 3 2 1; 3 4 1; 3 5 1; 4 2 1; 4 3 1; 4 6 1; 5 3 1; 5 6 1; 6 4 1; 6 5 1]); xy = [0 1; 1 2; 1 0; 2 0; 2 2; 3 1]; % 每一個列向量是一組 (x, y) 座標 gplot(A, xy, '-o') % 畫出無向圖（Undirected Graph）

22 Adjacency Matrix Real-world examples 台北捷運：n>100

23 稀疏矩陣有趣的例子-1 Bucky 球，此圖包含了 60 個三度空間中的點，每一點和他的三個鄰近點都是等距離，可用 bucky 指令來產生這些點的鄰近矩陣，並用 gplot 來顯示圖形 範例5-7：gplot02.m [A, xy] = bucky; % A 為鄰近矩陣，xy 為座標 gplot(A, xy, '-o'); % 畫出無向圖（Undirected Graph） axis equal % 設定 x 軸和 y 軸的刻度一樣）

24 畫出一棵電腦圖學中的樹 treeplot指令可畫出一棵電腦圖學中的樹 範例5-8：treePlot01.m
使用 parent 向量來代表這一棵樹，其中 parent(1)=0 則代表第一個節點是此樹的根節點（Root），而 parent(i)=j 代表第 i 個節點的父親是第 j 個節點，例如 parent(5)=4 代表第5個節點的父親是第 4 個節點，依此類推 parent = [ ]; treeplot(parent);

25 稀疏矩陣有趣的例子-2 是由 NASA （美國太空總署）所主導的計畫，其中包含計算流過機翼的氣流所造成的作用力，由於必須進行偏微分方程的數值運算，所以必須對二維空間進行三角化切割，其鄰近矩陣即為一個稀疏矩陣，您可在 MATLAB 下執行 airfoil 指令即可產生相關圖形，相關說明則記載於 airfoil.m，可經由 type airfoil.m 來檢視之。

26 5-4 稀疏矩陣的運算 MATLAB 針對完全矩陣設計的運算與函式，也都適用於稀疏矩陣，而且輸出也是大部分以稀疏矩陣的方式來表示
5-4 稀疏矩陣的運算 MATLAB 針對完全矩陣設計的運算與函式，也都適用於稀疏矩陣，而且輸出也是大部分以稀疏矩陣的方式來表示 若計算過程包含稀疏及完全矩陣，則計算結果的表示方式就依情況而變，其規則可見 MATLAB 線上輔助說明

27 稀疏矩陣的運算 稀疏矩陣的運算還包含下列幾種：
排列（Permutation）及重排（Reordering） 因子分解（Factorization） 線性聯立方程式的求解 固有值及奇異值的計算 這方面的運算，牽涉到很多數學理論，有興趣的讀者，可參考 MATLAB 的手冊，或在 MATLAB 指令視窗下輸入「help sparfun」以取得線上輔助說明

28 5-5 本章指令彙整 指令 功能 B = sparse(A) 將一個完全矩陣 A 轉換成稀疏矩陣 B
5-5 本章指令彙整 指令 功能 B = sparse(A) 將一個完全矩陣 A 轉換成稀疏矩陣 B B = sparse(i, j, s, m , n ) 直接產生一個稀疏矩陣 B，其中 ： i 是列索引 j 是行索引 s 是非零元素形成的向量 m 是 B 的列維度 n 是 B 的行維度 i、j、s 都是長度相同的向量 B = spdiags(D, p, m, n) 由矩陣 D 的對角線元素來建構一個稀疏矩陣 B，其中： D 的每一個直行代表矩陣的對角線向量 p 代表對角線的位置（0 代表主對角線，-1 代表向下位移一單位的次對角線，1 代表向上位移一單位的次對角線，依此類推） m 與 n 則分別代表矩陣的列維度與行維度 full(B) 以完全矩陣來顯示矩陣 B

29 5-5 本章指令彙整 spconvert(C) 將一個 m×3 的矩陣 C，轉換成稀疏矩陣，其中： 第一直行代表列索引 第二直行代表行索引
5-5 本章指令彙整 spconvert(C) 將一個 m×3 的矩陣 C，轉換成稀疏矩陣，其中： 第一直行代表列索引 第二直行代表行索引 第三直行則是非零的元素值 nnz(C) 傳回稀疏矩陣 C 的非零元素個數 nonzeros(C) 傳回稀疏矩陣 C 的所有非零元素形成的行向量 nzmax(C) 傳回稀疏矩陣 C 的目前可容納非零元素個數的最大值，當nnz > nzmax 時，MATLAB 會動態調增配置記憶體給 nzmax，以儲存新增的非零元素 spy(C) 觀看稀疏矩陣 C 的非零元素分佈情況 gplot(A, xy, '-o'); 畫出無向圖（Undirected Graph） treeplot(nodes) 畫出樹狀圖