資料結構——矩陣壓縮與壓縮矩陣的轉置與乘法計算
阿新 • • 發佈：2018-12-15
為什麼要對矩陣進行壓縮儲存呢？對於一個n*m的矩陣，我們一般會想到開一個n*m的二維陣列來儲存，這樣計算操作都很方便模擬，但當一個矩陣很大時，這樣對於空間的開銷與浪費是很可怕的，尤其是當矩陣變成多維時。但我們往往不會在矩陣每一個位置都存有資料，很多矩陣元素其實是0，我們需要記錄的只是那些非零元素，於是我們可以記錄非零元素的位置與值，這樣便可以大大減少空間上的浪費。
矩陣壓縮儲存程式碼（注意要按行順序遍歷儲存矩陣非零元素）：
1 typedef struct { 2 int i, j; //非零元行下標與列下標 3 ElemType e; 4 }Triple;5 6 typedef struct { 7 Triple data[MAXSIZE + 1]; 8 int mu, nu, tu; //行數， 列數， 非零元個數 9 }TSMatrix;
矩陣轉置是矩陣運算中極為重要的操作，但矩陣壓縮後我們不能再按原來的方式用n*m的遍歷方式進行轉置操作，那我們應該怎麼進行轉置呢，首先想到的是以列為參考值對元素進行遍歷轉置，這樣訪問元素的順序正好是轉制後儲存的順序，程式碼如下：
1 Status TransposeSMatrix(TSMatrix M, TSMatrix & T) { 2 T.mu = M.nu;3 T.nu = M.mu; 4 T.tu = M.tu; 5 if(T.tu) { 6 int q = 1; 7 for(int col = 1; col <= M.nu; ++col) 8 for(int p = 1; p <= M.tu; ++p) 9 if(M.data[p].j == col) { 10 T.data[q].i = M.data[p].j; 11 T.data[q].j = M.data[p].i;12 T.data[q].e = M.data[p].e; 13 ++q; 14 } 15 } 16 return OK; 17 }
我們注意到當矩陣不那麼稀疏時，比如如果這是一個稠密矩陣甚至n*m的空間裡每一個位置都是一個非零元素，這時這個程式碼的複雜度達到了O(mu*nu*nu)，比直接用陣列儲存的O(mu*nu)還要大一個數量級，這樣節省一點空間浪費大量時間的做法顯然是不合適的，我們考慮改進演算法，我們注意到如果我們預先知道矩陣M每一列的第一個非零元素在T中應有的位置，那麼轉置中我們就可以直接把他放在那個位置，這時，我們考慮開輔助陣列記錄每一列第一個非零元素的位置，以及每一列非零元素的數量，但其實，每一列的第一個非零元素的位置就是前一列第一個非零元素位置加上前一列非零元素數量，我們每次維護剩餘元素在剩餘矩陣中每一列第一個非零元素的位置陣列，就得到了下面的O(nu + tu)程式碼:
1 Status FastTransposeSMtrix(TSMatrix M, TSMatrix &T) { 2 T.mu = M.nu; 3 T.nu = M.mu; 4 T.tu = M.tu; 5 if(T.tu) { 6 int *num = new int[M.nu + 1], *cpot = new int[M.tu + 1]; 7 for(int col = 1; col <= M.nu; ++ col) 8 num[col] = 0; 9 for(int t = 1; t <= M.tu; ++t) ++num[M.data[t].j]; 10 cpot[1] = 1; 11 for(int col = 2; col <= M.nu; ++ col) 12 cpot[col] = cpot[col - 1] + num[col - 1]; 13 for(int p = 1; p <= M.tu; ++p) { 14 int col = M.data[p].j, q = cpot[col]; 15 T.data[q].i = M.data[p].j; 16 T.data[q].j = M.data[p].i; 17 T.data[q].e = M.data[p].e; 18 ++cpot[col]; 19 } 20 delete num; 21 delete cpot; 22 } 23 return OK; 24 }
矩陣轉置到此就結束了，寫進一個cpp效果就像下面這樣，感興趣的可以自己測試下效果：
1 #include <iostream> 2 #include <cstdio> 3 #include <cstdlib> 4 #include <cmath> 5 #include <algorithm> 6 #include <cstring> 7 #include <vector> 8 #include <string> 9 #include <queue> 10 #include <map> 11 #include <set> 12 13 #define FRER() freopen("in.txt", "r", stdin); 14 #define INF 0x3f3f3f3f 15 16 using namespace std; 17 18 //函式狀態碼定義 19 #define TRUE 1 20 #define FALSE 0 21 #define OK 1 22 #define ERROR 0 23 #define INFEASIBLE -1 24 //#define OVERFLOW -2 25 26 typedef int Status; 27 typedef int ElemType; 28 29 #define MAXSIZE 12500 30 31 typedef struct { 32 int i, j; //非零元行下標與列下標 33 ElemType e; 34 }Triple; 35 36 typedef struct { 37 Triple data[MAXSIZE + 1]; 38 int mu, nu, tu; //行數， 列數， 非零元個數 39 }TSMatrix; 40 41 Status TransposeSMatrix(TSMatrix M, TSMatrix & T) { 42 T.mu = M.nu; 43 T.nu = M.mu; 44 T.tu = M.tu; 45 if(T.tu) { 46 int q = 1; 47 for(int col = 1; col <= M.nu; ++col) 48 for(int p = 1; p <= M.tu; ++p) 49 if(M.data[p].j == col) { 50 T.data[q].i = M.data[p].j; 51 T.data[q].j = M.data[p].i; 52 T.data[q].e = M.data[p].e; 53 ++q; 54 } 55 } 56 return OK; 57 } 58 59 Status FastTransposeSMtrix(TSMatrix M, TSMatrix &T) { 60 T.mu = M.nu; 61 T.nu = M.mu; 62 T.tu = M.tu; 63 if(T.tu) { 64 int *num = new int[M.nu + 1], *cpot = new int[M.tu + 1]; 65 for(int col = 1; col <= M.nu; ++ col) 66 num[col] = 0; 67 for(int t = 1; t <= M.tu; ++t) ++num[M.data[t].j]; 68 cpot[1] = 1; 69 for(int col = 2; col <= M.nu; ++ col) 70 cpot[col] = cpot[col - 1] + num[col - 1]; 71 for(int p = 1; p <= M.tu; ++p) { 72 int col = M.data[p].j, q = cpot[col]; 73 T.data[q].i = M.data[p].j; 74 T.data[q].j = M.data[p].i; 75 T.data[q].e = M.data[p].e; 76 ++cpot[col]; 77 } 78 delete num; 79 delete cpot; 80 } 81 return OK; 82 } 83 84 int main() 85 { 86 //FRER() 87 88 TSMatrix MatrixA, MatrixB; 89 cin >> MatrixA.mu >> MatrixA.nu >> MatrixA.tu; 90 for(int i = 1; i <= MatrixA.tu; ++i) 91 cin >> MatrixA.data[i].i >> MatrixA.data[i].j >> MatrixA.data[i].e; 92 //TransposeSMatrix(MatrixA, MatrixB); 93 FastTransposeSMtrix(MatrixA, MatrixB); 94 for(int i = 1; i <= MatrixB.tu; ++i) 95 cout << MatrixB.data[i].i << ' ' << MatrixB.data[i].j << ' ' << MatrixB.data[i].e << endl; 96 97 return 0; 98 } 99 100 /*測試資料 101 6 7 8 102 1 2 12 103 1 3 9 104 3 1 -3 105 3 6 14 106 4 3 24 107 5 2 18 108 6 1 15 109 6 4 -7 110 */
接下來是壓縮矩陣的乘法，類似轉置，我們直接在矩陣結構體裡開一個數組用來存矩陣每一行第一個非零元素的在M中位置，就像這樣：
1 typedef struct { 2 Triple data[MAXSIZE + 1]; //非零元三元組表 3 int rops[MAXSIZE + 1]; //各行第一個非零元位置 4 int mu, nu, tu; //行數，列數，非零元個數 5 }RLSMatrix;
由於稀疏矩陣相乘不一定還是稀疏矩陣，所以我們要根據結果判斷元素是否是非零元，模擬乘法運算如下：
1 Status MultSMatrix(RLSMatrix M, RLSMatrix N, RLSMatrix & Q) { 2 if(M.nu != N.mu) return ERROR; 3 Q.mu = M.mu; 4 Q.nu = N.nu; 5 Q.tu = 0; 6 int *ctemp = new int[Q.nu + 1]; 7 8 if(M.mu * N.nu != 0) { 9 for(int arow = 1; arow <= M.mu; ++arow) { 10 int tp; 11 memset(ctemp, 0, (Q.nu + 1) * sizeof(int)); 12 Q.rops[arow] = Q.tu + 1; 13 if(arow < M.mu) 14 tp = M.rops[arow + 1]; 15 else 16 tp = M.tu + 1; 17 for(int p = M.rops[arow]; p < tp; ++p) { //對當前行的每一個非零元進行計算 18 int brow = M.data[p].j; //找到對應元在N中的標號 19 int t; 20 if(brow < N.mu) t = N.rops[brow + 1]; 21 else t = N.tu + 1; 22 for(int q = N.rops[brow]; q < t; ++q) { 23 int ccol = N.data[q].j; 24 ctemp[ccol] += M.data[p].e * N.data[q].e; 25 } 26 } 27 for(int ccol = 1; ccol <= Q.nu; ++ccol) //壓縮儲存該行非零元 28 if(ctemp[ccol]) { 29 if(++Q.tu > MAXSIZE) return ERROR; 30 Q.data[Q.tu].i = arow; 31 Q.data[Q.tu].j = ccol; 32 Q.data[Q.tu].e = ctemp[ccol]; 33 } 34 } 35 } 36 return OK; 37 }