动态规划法——矩阵连乘问题

给定n个矩阵｛A1,A2,...,An｝，其中Ai与Ai+1是可乘的，i=1，2...，n-1。如何确定计算矩阵连乘积的计算次序，使得依此次序计算矩阵连乘积需要的数乘次数最少。

按设计动态规划算法的步骤解题。
(1)找出最优解的性质，并刻划其结构特征。
(2)递归地定义最优值。
(3)以自底向上的方式计算出最优值。
(4)根据计算最优值时得到的信息，构造最优解(由子结构的最优解得到原先大问题的最优解)。

问题辅助分析：

（1）问题的描述

给定n个矩阵{A1,A2,…,An}，其中Ai与Ai+1是可乘的，i=1,2,…,n-1。要算出这n个矩阵的连乘积A1A2…An。由于矩阵乘法满足结合律，故计算矩阵的连乘积可以有许多不同的计算次序。这种计算次序可以用加括号的方式来确定。若一个矩阵连乘积的计算次序完全确定，也就是说该连乘积已完全加括号，则可以依此次序反复调用2个矩阵相乘的标准算法计算出矩阵连乘积。

完全加括号的矩阵连乘积可递归地定义为： （1）单个矩阵是完全加括号的； （2）矩阵连乘积A是完全加括号的，则A可表示为2个完全加括号的矩阵连乘积B和C的乘积并加括号，即A=(BC)。

 例如，矩阵连乘积A1A2A3A4有5种不同的完全加括号的方式：（A1（A2（A3A4））），（A1（（A2A3）A4）），（（A1A2）（A3A4）），（（A1（A2A3））A4），（（（A1A2）A3）A4）。每一种完全加括号的方式对应于一个矩阵连乘积的计算次序，这决定着作乘积所需要的计算量。若A是一个p×q矩阵，B是一个q×r矩阵，则计算其乘积C=AB的标准算法中，需要进行pqr次数乘。

 为了说明在计算矩阵连乘积时，加括号方式对整个计算量的影响，先考察3个矩阵{A1,A2,A3}连乘的情况。设这三个矩阵的维数分别为10×100，100×5，5×50。加括号的方式只有两种：（（A1A2）A3），（A1（A2A3）），第一种方式需要的数乘次数为10×100×5＋10×5×50＝7500，第二种方式需要的数乘次数为100×5×50＋10×100×50＝75000。第二种加括号方式的计算量时第一种方式计算量的10倍。由此可见，在计算矩阵连乘积时，加括号方式，即计算次序对计算量有很大的影响。于是，自然提出矩阵连乘积的最优计算次序问题，即对于给定的相继n个矩阵{A1,A2,…,An}（其中矩阵Ai的维数为pi-1×pi，i＝1,2,…,n），如何确定计算矩阵连乘积A1A2…An的计算次序（完全加括号方式），使得依此次序计算矩阵连乘积需要的数乘次数最少。

穷举搜索法的计算量太大，它不是一个有效的算法，本实验采用动态规划算法解矩阵连乘积的最优计算次序问题。

（2）算法设计思想

动态规划算法的基本思想是将待求解问题分成若干个子问题，先求解子问题，然后从这些子问题的解得到原问题的解。与分治法不同的是，动态规划法经分解得到的子问题往往不是相互独立的，前一子问题的解为后一子问题的解提供有用的信息，可以用一个表来记录所有已解决的子问题的答案，不管该子问题以后是否被用到，只要它被计算过，就将其结果填入表中。

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

#include<iostream>
#include<iomanip> 
#include<memory.h>
using namespace std;
//求出矩阵A[i:j]的最少数乘次数m[i][j]，和记录矩阵A[i:j]此时的断开位置s[i][j]
void MatrixChain(int *p, int n, int m[100][100], int s[100][100] )
{
        for (int i = 1; i <= n; i++) m[i][i] = 0; //初始化，使用矩阵的下标从1，1开始        
		for (int r = 2; r <= n; r++)
           for (int i = 1; i <= n - r+1; i++) {
              int j=i+r-1;
              m[i][j] = m[i+1][j]+ p[i-1]*p[i]*p[j];
              s[i][j] = 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;}
              }
          }
}

//利用断开位置s[i][j]输出矩阵A[i:j]的最优加括号方式
void print_optimal(int s[100][100], int i ,int j, int a[])
{
	if(i==j)	cout<<" A["<<a[i-1]<<","<<a[i]<<"]";
	else{
		cout<<" ( ";
		print_optimal(s,i,s[i][j],a);
		print_optimal(s,s[i][j]+1,j,a);
		cout<<" ) ";
	}
}

void putout(int a[100][100], int n)
{
	int i,j;
	for(i=1;i<=n;i++){
		for(j=1;j<=n;j++)
			cout<<setw(6)<<a[i][j];
		cout<<endl;
	}
}

int main()  
{   
	int n,k; 
	int array[100], m[100][100], s[100][100];
	memset(array,0,sizeof(array));
	memset(m,0,sizeof(m));
	memset(s,0,sizeof(s));

	while(1)
	{
		cout<<endl<<"输入矩阵的个数 "<<endl; 
		cin>>n;  
		if(n<=0)break;
		
		cout<<"请输入 "<<n+1<<" 个正整数，分别是各个矩阵的行列数"<<endl;  
		for(k=0;k<n+1;k++)  
			cin>>array[k];  
		cout<<endl; 
		cout<<"原始数据为以下矩阵"<<endl;
		for(k=0;k<n;k++)
			cout<<endl<<setw(3)<<array[k]
			<<" * "<<setw(3)<<array[k+1];  
		cout<<endl; 

		MatrixChain(array,n,m,s); 
		
		cout<<"m矩阵"<<endl; 
		putout(m,n);

		cout<<"s矩阵"<<endl; 
		putout(s,n);
		
		cout<<endl<<"最优的运算方式的乘法次数为: "<<m[1][n]<<endl
			<<"加括号的方式为："<<endl;
		print_optimal(s,1,n,array);
		
	}
	return 0;
}