程序算法藝術與實踐:遞歸策略之矩陣乘法問題 · 程序算法藝術與實踐

**矩陣預備知識** 在數學中，矩陣（Matrix）是指縱橫排列的二維數據表格，最早來自于方程組的系數及常數所構成的方陣。這一概念由19世紀英國數學家凱利首先提出。矩陣是高等代數學中的常見工具，也常見于統計分析等應用數學學科中。并且在ACM競賽，有很多涉及到矩陣知識的題。許多算法都會結合矩陣來處理，而比較具有代表性的矩陣算法有：矩陣快速冪、高斯消元等等。例如下面的圖片就是一個矩陣： ![](https://box.kancloud.cn/e94ce4faf3a86cf2cbf59ca969c29348_131x82.jpg) 上述矩陣是一個 4 × 3 矩陣：某矩陣 A 的第 i 行第 j 列，或 I , j位，通常記為 A[i,j] 或 Aij。在上述例子中A[3,3]=37。此外 Aij = (aij)，意為 A[i,j] = aij。①? 矩陣相乘的規則：矩陣與矩陣相乘第一個矩陣的列數必須等于第二個矩陣的行數假如第一個是m*n的矩陣第二個是n*p的矩陣則結果就是m*p的矩陣且得出來的矩陣中元素具有以下特點：第一行第一列元素為第一個矩陣的第一行的每個元素和第二個矩陣的第一列的每個元素乘積的和以此類推第i行第j列的元素就是第一個矩陣的第i行的每個元素與第二個矩陣第j列的每個元素的乘積的和。②? 單位矩陣： n*n的矩陣? mat ( i , i )=1;? 任何一個矩陣乘以單位矩陣就是它本身 n*單位矩陣=n，可以把單位矩陣等價為整數1。（單位矩陣用在矩陣快速冪中）例如下圖就是一個3*3的單位矩陣： ![](https://box.kancloud.cn/1bcdd3ce798d1244d462b1dc3dba3e35_84x66.jpg) ### 問題描述: 設A1,A2,…,An為矩陣序列，Ai為Pi-1×Pi階矩陣，i = 1,2,…,n. 確定乘法順序使得元素相乘的總次數最少.輸入：向量P = <P0, P1, … , Pn> **實例：?** P = <10, 100, 5, 50> ?A1: 10 × 100, A2: 100 × 5, A3: 5 × 50 **乘法次序：** (A1 A2)A3: 10 × 100 × 5 + 10 ×5 × 50 = 7500 ? ? ? A1(A2 A3): 10 × 100 × 50 + 100 × 5 × 50 = 75000 搜索空間的規模先將矩陣鏈加括號分為兩部分，即P=A1*A2*...*An=（A1*A2...*Ak）*（Ak+1*...An），則有f(n)=f(1)*f(n-1)+f(2)*f(n-2)+...+f(n-1)*f(1)種方法。f(n)為一個[Catalan數](http://baike.baidu.com/view/1163998.htm)，所以一般的方法要計算![](https://box.kancloud.cn/fdde16a514a90a7476949b6e5a2018e3_59x40.png) 種。 ### 動態規劃算法輸入P=< P0, P1, …, Pn>，Ai..j 表示乘積 AiAi+1…Aj 的結果，其最后一次相乘是： ![](https://box.kancloud.cn/5e3f8782b424279ccac68470d60cd247_123x20.png) m[i,j] 表示得到Ai..j的最少的相乘次數。遞推方程： ![](https://box.kancloud.cn/eb9cf4224c1dcf719446d6d31757c1d0_338x55.png) 為了確定加括號的次序，設計表s[i,j],記錄求得最優時最一位置。 ### 算法遞歸實現由上面的遞歸公式，很容易得到算法的遞歸實現,用一個N=10, P=<30,35,15,5,40,20,10,8,9,60,80>的簡單實例.參考代碼如下所示： ~~~ int RecurMatrixChain(int *ptrArray,int start,int end){ m[start][end]=100000; s[start][end]=start; if(start==end) m[start][end]=0; else{ for(int k=start;k<end;k++){ int q=RecurMatrixChain(ptrArray,start,k)+RecurMatrixChain(ptrArray,k+1,end)+ptrArray[start]*ptrArray[k+1]*ptrArray[end+1]; if(q<m[start][end]){ m[start][end]=q; s[start][end]=k; } } } return m[start][end]; } ~~~ 測試結果如下所示: ![](https://box.kancloud.cn/e352c4594c636f5126721b12a3e83095_640x113.jpg) ### 遞歸實現的復雜性復雜性滿足遞推關系： ![](https://box.kancloud.cn/caa43bbc9450b965dee934cdf2e9823e_346x115.png) 由數學歸納法可得：![](https://box.kancloud.cn/157b29efc88299b5886fda1bfd198b4a_74x23.png) 可見遞歸實現的復雜性雖然較一般算法有改進，但還是較高。分析原因，主要是**子問題重復程度高**。如下圖所示(ps:此圖來自網絡)： ![](https://box.kancloud.cn/a53b11f406a37412100e1429fd1246ad_679x269.png) 1..4表示計算Ai..j中i=1,j=4的子問題，其子問題包括A1..1，而A1..2，A1..3中都包括子問題A1..1，所以很多子問題被重復計算了多次。于是，我們想到用自底向上的迭代實現。 ### 算法迭代實現迭代實現主要思想是子問題由小到大，每個子問題只計算一次，并且把結果保存起來，后來用到這個子問題時，直接代入。 ~~~ void MatrixChain(int *P,int n){ int r,i,j,k,t; for(i=0;i<N;i++) for(j=0;j<N;j++) m[i][j]=0; //r為當前計算的鏈長（子問題規模） for(r=2;r<=n;r++){ //n-r+1為最后一個r鏈的前邊界 for(i=0;i<n-r+1;i++){ //計算前邊界為r，鏈長為r的鏈的后邊界 j=i+r-1; //將鏈ij劃分為A(i) * ( (A(i+1) ... A(j) ) m[i][j]=m[i+1][j]+P[i]*P[i+1]*P[j+1]; //記錄分割位置 s[i][j]=i; for( k=i+1;k<j-1;k++){ //將鏈ij劃分為( A(i)...A(k) )* ( (A(k+1) ... A(j) ) t=m[i][k]+m[k+1][j]+P[i]*P[i+1]*P[j+1]; if(t<m[i][j]){ m[i][j]=t; s[i][j]=k; } } } } } ~~~ **其中: 迭代實現的復雜性.**行7,9,16的循環為O(n)，外層循環為O(1)，所以算法復雜度W(n)=O(n^3).**迭代過程的一個實例**:子問題由小到大的計算過程如下圖所示(ps:此圖來自網上)： ![](https://box.kancloud.cn/2e909f17535989b94365d004cf034f59_583x379.jpg) 再寫一個打印結果，以及打印優化函數備忘錄m和標記函數的s的函數,： ~~~ void PrintMatrixChain(int s[][N],int i,int j){ if (i==j) { cout<<"A"<<i+1; } else { cout<<"("; PrintMatrixChain(s, i, s[i][j]); PrintMatrixChain(s, s[i][j]+1, j); cout<<")"; } } void PrintMS(int m[][N],int s[][N],int N){ for(int r=0;r<N;r++){ for(int i=0;i<N-r;i++){ int j=i+r; cout<<"m["<<i+1<<","<<j+1<<"]="<<m[i][j]<<"\t"; } cout<<endl; } for(int r=1;r<5;r++){ for(int i=0;i<N-r;i++){ int j=i+r; cout<<"s["<<i+1<<","<<j+1<<"]="<<s[i][j]+1<<"\t"; } cout<<endl; } } ~~~ **一個簡單的測試實例**:用一個N=6，P=<{8,30,80,35,5,10,20}>的簡單實例，運行上述代碼： ![](https://box.kancloud.cn/f8f041378b926670209bf80efc2f0a2b_641x388.jpg) ### HDOJ 4920 Matrix multiplication **Problem Description** Given two matrices A and B of size n×n, find the product of them.bobo hates big integers. So you are only asked to find the result modulo 3. **Input** The input consists of several tests. For each tests: The first line contains n (1≤n≤800). Each of the following n lines contain n integers — the description of the matrix A. The j-th integer in the i-th line equals Aij. The next n lines describe the matrix B in similar format (0≤Aij,Bij≤109). **Output** For each tests: Print n lines. Each of them contain n integers — the matrix A×B in similar format. **Sample Input** 1 0 1 2 0 1 2 3 4 5 6 7 **Sample Output** 0 0 1 2 1 題目大意：求兩個矩陣相乘mod3的結果。這道題就是一個赤裸裸的矩陣相乘的題。但是要注意題目的時間復雜度，一定要進行優化。并且，此題還有二個坑點，如果把握不好就會超時。① 就是Mod 3 時，一定不能在矩陣相乘算法的for循環里，如果放在相乘算法里會超時。②在進行相乘之前把所有元素先Mod 3 這樣做也會優化一定的時間。因為題目所給數據并不是很大，所以即使把mod 3 放到結束輸出語句里面就可以了，不用擔心相乘的過程會超出int型。 **源代碼參考**： ~~~ int X[MAXN][MAXN]; int Y[MAXN][MAXN]; int t[MAXN][MAXN]; int main() { while(scanf("%d",&N)!=EOF){ for(int i=0;i<N;i++){ for(int j=0;j<N;j++){ scanf("%d",&X[i][j]); X[i][j]%=3; } } for(int i=0;i<N;i++){ for(int j=0;j<N;j++){ scanf("%d",&Y[i][j]); Y[i][j]%=3; t[i][j]=0; } } for(int i=0;i<N;i++){ for(int k=0;k<N;k++) if(X[i][k]) for(intj=0;j<N;j++) t[i][j]=(t[i][j]+X[i][k]*Y[k][j]%3)%3; } for(int i=0;i<N;i++){ printf("%d",t[i][0]); for(int j=1;j<N;j++){ printf("%d",t[i][j]); } printf("\n"); } } return 0; } ~~~ 單獨考察矩陣相乘的題目并不多見，普遍都是結合矩陣快速冪，高斯消元等矩陣算法。進行這些算法之前我們必須要掌握矩陣相乘算法.我們需掌握簡單的相乘規則，才能學習后面的一些矩陣算法。同時，為了以后學習算法以及做題的需要，我們也要記得矩陣相乘時的算法復雜度及優化細節。關于[程序算法藝術與實踐](http://blog.csdn.net/column/details/tac-programalgrithm.html)更多討論與交流，敬請關注本博客和新浪微博[songzi_tea](http://weibo.com/songzitea).