搞定矩阵连乘的算法技巧

        给定n个矩阵，其中Ai与Ai+1是可乘的（i=1，2，…… n-1）。考察这n个矩阵的连乘积A1A2A2……An。如何确定计算矩阵连乘积的计算次序，使得依此次序计算矩阵连乘积需要的数乘次数最少。

        由于矩阵乘法满足结合律，所以计算矩阵的连乘可以有许多不同的计算次序。这种计算次序可以用加括号的方式来确定。 若一个矩阵连乘积的计算次序完全确定，也就是说该连乘积已完全加括号，则可以依此次序反复调用2个矩阵相乘的标准算法计算出矩阵连乘积。完全加括号的矩阵连乘积可递归的定义为：

        1、单个矩阵是完全加括号的；

        2、矩阵连乘积A是完全加括号的，则A可表示为两个完全加括号的矩阵连乘积B和C的乘积并加括号，即A=（BC）。

         例如，矩阵连乘积可以有以下5种完全加括号方式：

,,,,

         矩阵A和B可乘的条件是矩阵A的列数等于矩阵B的行数，若A是一个p*q矩阵，B是一个q*r矩阵，则其乘积C=AB是一个p*r的矩阵，共需要p*q*r次数乘。

        穷举法：列举出所有可能的计算次序，并计算出每一种计算次序相应需要的数乘次数，从中找出一种数乘次数最少的计算次序。 

/*穷举法：列举出所有可能的计算次序，并计算出每一种计算次序相应需要的数乘次数，从中找出一种数乘次数最少的计算次序*/
void matrixMultiply(int** a, int** b, int** c, int ra, int ca, int rb, int cb)
/*
  ra、ca表示矩阵A的行数和列数，rb、cb表示矩阵B的行数和列数
*/
{
	if (ca != cb)
		//矩阵A和B可乘的条件是矩阵A的列数等于矩阵B的行数
		cout<<"矩阵不可乘"<<endl;
	//矩阵A的每一行分别与矩阵B的每一列相乘
	for (int i = 0; i < ra; i++)
	{
		for (int j = 0; j < cb; j++)
		{
			int sum = a[i][0] * b[0][j];
			//将A的每一行与B的每一列对应位置相乘并相加
			for (int k = 1; k < ca; k++)
				sum += a[i][k] * b[k][j];
			//C=A*B，记录乘积结果
			c[i][j] = sum;
		}
	}
}

算法复杂度分析：

       对于n个矩阵的连乘积，设其不同的计算次序为P(n)。 由于每种加括号方式都可以分解为两个子矩阵的加括号问题：(A1...Ak)(Ak+1…An)可以得到关于P(n)的递推式如下： 

 动态规划算法：

        将矩阵连乘积   简记为A[i:j] ，这里i≤j。考察计算A[i:j]的最优计算次序。设这个计算次序在矩阵 Ak和Ak+1之间将矩阵链断开，i≤k<j,则其相应完全加括号方式

计算量：A[i:k]的计算量加上A[k+1:j]的计算量，再加上 A[i:k]和A[k+1:j]相乘的计算量。

步骤：

       1、分析最优解的结构

        特征：计算A[i:j]的最优次序所包含的计算矩阵子链 A[i:k]和A[k+1:j]的次序也是最优的。 矩阵连乘计算次序问题的最优解包含着其子问题的最优解。这种性质称为最优子结构性质。问题的最优子结构性质是该问题可用动态规划算法求解的显著特征。

         2、建立递归关系

         设计算A[i:j]，1≤i≤j≤n，所需要的最少数乘次数m[i,j]，则原问题的最优值为m[1,n]

        当i=j时，A[i:j]=Ai，因此，m[i,i]=0，i=1,2,…,n

        当i<j时，

        3、计算最优值

        用动态规划算法解此问题，可依据其递归式以自底向上的方式进行计算。在计算过程中，保存已解决的子问题答案。每个子问题只计算一次，而在后面需要时只要简单查一下，从而避免大量的重复计算，最终得到多项式时间的算法。

void MatrixChain(int* p, int n, int** m, int** s)
{
	for (int i = 1; i <= n; i++)
		m[i][i] = 0;//自己与自己的相乘次数是0
	for (int r = 2; r <= n; r++)
		/* r表示当前相乘的矩阵个数 */
	{
		for (int i = 1; i <= n - r + 1; i++)
			/* i表示当前相乘的起始矩阵，j表示当前相乘的终止矩阵*/
		{
			int j = i + r - 1;
			m[i][j] = m[i][i] + m[i + 1][j] + p[i - 1] * p[i] * p[j];
			s[i][j] = i;//在i之后断开
			for (int k = i + 1; k < j; k++)
			{
				int t = m[i][k] + m[k + 1][j] + p[i - 1] * p[k] * p[j];//找下一个断开位置
				if (t < m[i][j])
				{
					m[i][j] = t;//记录最优值
					s[i][j] = k;//在k之后断开
				}
			}
		}
	}
}

            4、构造最优解

void Traceback(int i, int j, int** s)
{
	if (i == j)
		return;
	Traceback(i, s[i][j], s);
	Traceback(s[i][j] + 1, j, s);
	cout << "Multiply A" << i << "," << s[i][j];
	cout << "and A " << (s[i][j] + 1) << "," << j << endl;
}

           s[i][j]中的数表明，计算矩阵链A[i:j]的最佳方式是在矩阵A(k)和A(k+1)之间断开，即最优加括号方式为。