矩阵连乘
矩阵连乘详解html
--crystal yi数组
既然这篇文章叫作矩阵连乘详解，那么我就不能辜负详解这两个字，只有把一个原来不懂的的人弄懂了，才叫详解。函数
言归正传，首先让咱们复习一下矩阵连乘的有关知识。对于矩阵知识很了解的人能够跳过矩阵知识这块内容，不过笔者建议最好复习一下。大数据
矩阵知识:url
矩阵的乘法：spa
左面的矩阵的行数决定结果的行数，右面的矩阵的列数决定结果的行数，两个矩阵在相乘时必定要有相同的行数或者列数。code
矩阵相乘只有在第一个矩阵的列数和第二个矩阵的行数相同時才有定义。假如A为m×n矩阵，B为n×p矩阵，则他們的乘AB(有时记作A·B)会是一个m×p矩阵。xml
例如一个2x3的矩阵与一个3x2的矩阵的乘会是一个2x2的矩阵 。htm
例子：blog
矩阵连乘：
设有矩阵M1,M2,M3,M4，
其维数分别是10×20, 20×50, 50×1 和1×100，现要求出这4个矩阵相乘的结果。咱们知道，若矩阵A的维数是p×q，矩阵B的维数是q×r，则A与B相乘后所得矩阵AB的维数是p×r。按照矩阵相乘的定义，求出矩阵AB中的一个元素须要作q次乘法（及q-1次加法）。这样，要计算出AB就须要作p×q×r次乘法。为简单起见，且因为加法比一样数量的乘法所用时间要少得多，故这里咱们暂不考虑加法的计算量。因为矩阵连乘知足结合律，故计算矩阵连乘的方式能够有多种。
例如，咱们能够按M1(M2(M3M4))的方式去计算，
也能够按(M1(M2M3))M4的方式去计算，所得结果是相同的。
可是值得注意的是，
按前一方式计算须要作125,000次乘法，
而按后一方式计算只须要作2,200次乘法。
因而可知，矩阵连乘的运算次序对于所须要的计算量
（所需乘法次数）有着极大的影响。
M3M4：50*1*100=5,000；M2(M3M4)：20*50*100=100,000
M1(M2(M3M4))：10*20*100=20,000
(M2M3)：20*50*1=1000；(M1(M2M3))：10*20*1=200 ；
(M1(M2M3))M4：10*1*100=1000
如何解决矩阵连乘问题：
分析：
设要求出矩阵连乘MiMi+1……Mj-1Mj（i<j）所需的最少乘法次数。
因共有j-i+1个矩阵，故称这个矩阵连乘的规模是j-i+1
按照作最后一次乘法的位置进行划分，该矩阵连乘一共可分为j-i种状况即有(j-i)种断开方式：Mi(Mi+1┅Mj)，(MiMi+1)(Mi+2┅Mj)，┅，(MiMi+1┅Mj-1)Mj。其中任一对括号内的矩阵个数（即规模）不超过j-i。若咱们已知任一个规模不超过j-i的矩阵连乘所需的最少乘法次数，咱们就能够很容易地计算出矩阵连乘MiMi+1┅Mj-1Mj（i<j）所需的最少乘法次数，其方法以下。将上述的j-i种状况表示为通式：(Mi┅Mk) (Mk+1┅Mj)（i≤k<j）。
记第t个矩阵Mt的列数为rt，并令rt-1为矩阵Mt的行数。
则Mi┅Mk连乘所得是ri-1×rk维矩阵，
Mk+1┅Mj连乘所得是rk×rj维矩阵，
故这两个矩阵相乘须要作ri-1×rk×rj次乘法
因为在此以前咱们已知
任一个规模不超过j-i的矩阵连乘所需的最少乘法次数，故(Mi┅Mk)和(Mk+1┅Mj)所需的最少乘法次数已知，将它们分别记之为mi,k和mk+1,j。
形为(Mi┅Mk)(Mk+1┅Mj)的矩阵连乘所需的最少乘法次数为：
mi,k + mk+1,j + ri-1×rk×rj。
对知足i≤k<j 的共j-i种状况逐一进行比较，咱们就能够获得
矩阵连乘MiMi+1┅Mj-1Mj（i<j）所需的最少乘法次数mi,j为：
mi,j=min { mi,k + mk+1,j + ri-1×rk×rj } （i≤k<j）
因而在初始时咱们定义mi,i=0（至关于单个矩阵的状况）
求m1,2：即i=1, j=2，
就是2个矩阵，无需划分 (k=1,由于i<=k<j)
m1,2=min{ m1,1 + m2,2 + ri-1×ri×ri+1
=ri-1×ri×ri+1=r0×r1×r2= 10×20×50=10000，
求m23：即i=2, j=3，故ri-1×ri×ri+1=r1×r2×r3=20×50×1=1000，
求m13：即i=1, j=3，min（i≤k<j）{mik+mk+1,j+ri-1×rk×rj}=
min{ (m11+m23+r0×r1×r3)(k=1)，(m12+m33+r0×r2×r3)(k=2)}=
min{(0+1000+200), (10000+0+500)}= min{1200,10500}=1200
矩阵在C语言中的表示方法：
A[i][j] (i为行，j为列)
拥有了上面的准备知识看下面的内容就会轻松不少。
注：下文所用的mi,j或者m[i][j]表示MiMi+1┅Mj-1Mj（i<j）所需的最少乘法次数
首先，让咱们本身思考一下怎么解决一个n个矩阵连乘的问题，也许有了上面的知识你会想到下面的递归的方法：
假设已经算好了Mi┅Mk和Mk+1┅Mj这两个子矩阵各自相乘的最小的次数，剩下的只须要算最后一次乘积的次数就能够了，也就是这里的两个子矩阵相乘所涉及的乘法次数ri-1×rk×rj而后加上前面所得的两个子矩阵各自所需的最小的相乘次数就能够了。这样就获得了上面所说的这个式子：mi,j = mi,k + mk+1,j + ri-1×rk×rj 而后依次判断在k等于某个数的时候所得的mi,j的值，剩下的只须要对全部的mi,j的值取一个最小值就ok了，因而咱们就获得了上面所说的这个式子： mi,j=min { mi,k + mk+1,j + ri-1×rk×rj } （i≤k<j）
然而这样作并不能求出mi,j的最小值，由于咱们不知道 mi,k 和 mk+1,j的值。可是，这并不会难倒咱们，由于咱们知道对于 mi,k 和 mk+1,j 的求解能够基于递归思想。方法和求解mi,j的最小值同样。这样咱们就获得了递归的式子：
下面代码的p[i]表示的是第i个矩阵的列数，固然他的行数就用i-1表示，呵呵。
1 下面代码的p[i]表示的是第i个矩阵的列数，固然他的行数就用i-1表示，呵呵。 2 3 int min(int a,intb){ 4 5 return a<b?a:b; 6 7 } 8 9 int d(int i,int j){ 10 11 int i,j,k,t,u; 12 13 if(i==j) return 0; //i==j时就是mi,i，也就是单个矩阵，固然，单个矩阵不存在乘积因此为0 14 15 if(i==j-1) return p[i-1]*p[i]*p[j]; //明显，这是相邻的两个矩阵的乘积 16 17 u=d(i,i)+d(i+1,j)+p[i-1]*p[i]*p[j]; //1 18 19 for(k=i+1;k<j;k++){ //2 20 21 t=d(i,k)+d(k+1,j)+p[i-1]*p[k]*p[j]; //2 22 23 u=min(t,u); //3 24 25 } 26 27 }
1 对于从矩阵i到矩阵j的连乘，当k==i时所得的mi,j
2 对于从矩阵i到矩阵j的连乘，当i<k<j时所得的 mi,j
注：1 ， 2 所说的结合起来就是 3 也就是上面蓝字部分的分析
3 就是mi,j=min { mi,k + mk+1,j + ri-1×rk×rj } （i≤k<j）
但是，有兴趣的话你能够算一下用递归作的效率，你会发现，复杂度是O(2^n)，这个复杂度太大了，稍大一点的数据就爆了不只timelimitexceeded并且也会memorylimitexceeded，因此咱们就必须换一种方法改用dp(动态规划)，可是，咱们的思路仍是同样的，换句话说，只不过把递归的形式改成dp的形式而已，没什么新的东西，若是说有的话也只是记录了中间状态。请看下面的代码：
1 int dp_matrix_multiply(){ 2 3 int i,j,r,k,tem; 4 5 n-=1; 6 7 for(i=1;i<=n;i++) //1 8 9 m[i][i]=0; //1 10 11 for(r=2;r<=n;r++){ //2 12 13 for(i=1;i<=n-r+1;i++){ //3 14 15 j=i+r-1; //3 16 17 m[i][j]=m[i+1][j]+p[i-1]*p[i]*p[j]; //4 18 19 for(k=i+1;k<j;k++){ //4 20 21 tem=m[i][k]+m[k+1][j]+p[i-1]*p[k]*p[j]; //4 22 23 m[i][j]=tem<m[i][j]?tem:m[i][j]; //4 24 25 } 26 27 } 28 29 } 30 31 return m[1][n]; 32 33 }
初看，是否是很头痛？？我刚开始接触的时候看到别人这么写也很头大，可是思考了近3个小时后终于理解了。
首先咱们看一下注释1，这部应该很好理解吧，就是初始化单个矩阵连乘的最小次数，固然单个矩阵不可能与其自身相乘因此赋值为0
而后是注释2，这个我想应该是初学者最难以理解的地方之一了吧(包括当时的我也同样,呵呵)，其实这里的r就是 the matrix train length，连乘的矩阵个数。这里从2开始，有两个缘由：首先单个矩阵的乘积(就是其自身，这一步已经在1中阐明了)已经初始化好了，也就是说咱们已经求得了r==1时的若干单个矩阵最小的乘积次数，接下去固然要求r==2时，若干个由两个矩阵组成的的子矩阵乘积的最小的次数，而后依次类推。其次，对于这里的r的解释咱们能够看一下这个状态方程： mi,j=min { mi,k + mk+1,j + ri-1×rk×rj } （i≤k<j）
若是咱们要知道mi,j的值，咱们就必须知道 mi,k 和 mk+1,j的值，而mi,k和mk+1,j的求解又是经过更小的子模块来实现，以此类推最后在mi,i，也就是单个矩阵的地方终止子问题的进一步求解，由于没有比单个矩阵更小的子问题了，讲到这里咱们不妨抽点时间再看一下上面讲的递推的方法，你会发现递归的实质就是求解更小的子模块，这不是和这里的dp求解子模块的方法同样吗。并且你会发现，注释4的部分和递归求解的2，3部分几乎如出一辙，那是否是一样的意思呢，相信本身吧，你理解的是正确的，可是，你要注意一下dp的方法优越于递归的方法的主要缘由就是将递归中的t=d(i,k)+d(k+1,j)+p[i-1]*p[k]*p[j];改为了tem=m[i][k]+m[k+1][j]+p[i-1]*p[k]*p[j];这样用数组记录状态节约了大量的时间，由于递归方法中d(i,k)这个函数老是被屡次调用，不能用递归的方法解决大数据的缘由就在于此。
理解了上面的内容注释3的意思也不难看懂了，i和j就是MiMi+1┅Mj-1Mj（i<j）中的i和j，至于i的值为何从1到n-r+1以及j的值为何是i+r-1不用我说你也应该理解，i从1到n-r+1就是r个矩阵连乘的第一个矩阵的下标(注：这里矩阵的下标从1开始),既然这样那么j就是r个矩阵相乘最后那个矩阵的下标。
呵呵，经过上面的详解应该理解了吧。
打铁趁热赶快作一下下面的题目：
http://poj.org/problem?id=1651 pku 1651
1 #include <stdio.h> 2 3 #include <string.h> 4 5 #include <stdlib.h> 6 7 #define len 105 8 9 int m[len][len],p[len]; 10 11 int n; 12 13 bool input(){ 14 15 int i; 16 17 if(scanf("%d",&n)==1&&n>0){ 18 19 for(i=0;i<n;i++) 20 21 scanf("%d",&p[i]); 22 23 return true; 24 25 } 26 27 else return 0; 28 29 } 30 31 int dp_matrix_multiply(){ 32 33 int i,j,r,k,tem; 34 35 n-=1; 36 37 for(i=1;i<=n;i++) 38 39 m[i][i]=0; 40 41 for(r=2;r<=n;r++){ 42 43 for(i=1;i<=n-r+1;i++){ 44 45 j=i+r-1; 46 47 //i==k 48 49 m[i][j]=m[i+1][j]+p[i-1]*p[i]*p[j]; 50 51 for(k=i+1;k<j;k++){ 52 53 tem=m[i][k]+m[k+1][j]+p[i-1]*p[k]*p[j]; 54 55 m[i][j]=tem<m[i][j]?tem:m[i][j]; 56 57 } 58 59 } 60 61 } 62 63 return m[1][n]; 64 65 } 66 67 int main(void){ 68 69 int ans; 70 71 while(input()){ 72 73 ans=dp_matrix_multiply(); 74 75 printf("%d\n",ans); 76 77 } 78 79 return 0; 80 81 }