动态规划之矩阵连乘

两个矩阵相乘：设A为ra*ca的矩阵，B为rb*cb的矩阵，C=A*B，求C

 

首先数组a[ra][ca]、b[rb][cb]、c[r][c]分别存储矩阵A，B，C(r=ra,cb=c)

if(ca!=rb)   两个矩阵不能相乘，

else  A的第一行与B的第一列对应位置相乘再相加

void MatrixMultiply(int **a,int **b,int **c,int ra,int ca,int rb, int cb){  //矩阵相乘
    if(ra != cb)
        cout<<"can not multiply!"<<endl;
    for(int i=0; i<ra; i++){
        for(int j=0; j<cb; j++){
           int sum = a[i][0]*b[0][j];
           for(int k=1; k<ca; k++){
                sum += a[i][k]*b[k][j];
           }
           c[i][j] = sum;
        }
    }
}

若A为pxq的矩阵，B为qxr的矩阵，则AB为pxr的矩阵，需要进行pxqxr次数乘。

多个矩阵连乘：给n个矩阵{A1,A2,...,An}，其中相邻的矩阵是可乘的，计算A1A2...An，因为矩阵乘法满足结合律，通过加括号改变计算次序可以减小其计算量，例如：A1、A2、A3分别为10*100、100*5和5*50的矩阵（（A1A2）A3）需要10*100*5+10*5*50=7500次数乘，（A1(A2A3)）需要100*5*50+10*100*50=75000次数乘。

问题：利用完全加括号的方式确定最优计算次序，使得数乘次数最少，并输出加括号后的结果。

输入：第一行输入参与乘积的矩阵的个数n；第二行依次输入矩阵的行、列（例：A1:30*35 A2:35*15 A3:15*5 A4:5*10 A5:10*20 A6:20*25<30 35 15 5 10 20 25>）用p[n]数组保存输入的值。

可以发现：矩阵个数+1=第二行输入的数的个数

穷举法计算量太大。用动态规划算法解矩阵连乘积的最优计算次序。步骤如下：

1. 分析最优解结构

动态规划第一步：刻画问题的最优解结构特征。

将AiAi+1...Aj简记为A[i:j]，即考察A[1:n]的最优计算次序。设从Ak和Ak+1之间断开结果最优，即加括号((A1...Ak)(Ak+1...An))，则所需数乘次数为A[1:k]的计算量+A[k+1:n]的计算量+A[1:k]结果*A[k+1:n]的计算量。该问题具有最优子结构性质：（即最优解包含着其子问题的最优解,但这个子问题是不独立的（区别于分治法））A[1:n]的最优次序包含其子链A[1:k]和A[k+1:n]的次序也是最优的。断开的位置k由数组s[i][j]保存

2. 建立递归关系

动态规划第二步：递归的定义最优值。

用数组m[i][j]保存最优数乘次数，则原问题的最优解为m[1][n]。

当i=j时，A[i][j]表示一个矩阵，不需要进行乘法运算，因此m[i][i]=0（i=1,2,3...,n）

当i<j时，当从k断开时结果最优，m[i][j]=min{ m[i][k]+m[k+1][j]+p[i-1]*p[k]*p[j] }<i<=k<j>

 

自底向上思想！！！：先找出两个矩阵相乘的最优值，需要遍历每两个相邻矩阵，再利用遍历过得相邻两个相乘的结果计算三个相乘的最优，以此类推，最终得到n个相乘的最优值 

3. 计算最优值

计算m[1:n]时，为避免子问题被重复计算多次（动态规划的又一个特征），在计算过程中保存已解决的子问题的答案，每个子问题计算一次，后面需要时只需简单的查看使用，从而避免重复计算。具体实现如下：

​
void MatrixChain(int *p,int n,int m[][N],int s[][N]) { //构造最优解，p存储数据的维数
    for(int i=1; i<=n; i++)
        m[i][i] = 0;//i==j时,p[i:j]为单一矩阵，数乘次数为0
    for(int r=2; r<=n; r++) {//从2个矩阵相乘到n个矩阵相乘，自底向上的思想
        for(int i=1; i<=n-r+1; i++) { 
//若r=2即有两个矩阵相乘，将每两个相邻的矩阵相乘的最优值存入m数组
//当r=3时，可以利用已存入m的值得到不同的3个矩阵连乘的最优值，依次类推
            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;
                }
            }
        }
    }
}
//r=2时，其中A1A2计算次数存在m[1][2]中，对应位置s[1][2]=1;A2A3存在m[2][3]中，s[2][3]=2...
​

当r=n时，将得到最终的最优解

 4. 构造最优解

动态规划第四步：构造问题最优解，通过上面的函数，最优值被计算出来，也记录了矩阵链断开的位置，接下来就是要根据数组s中记录的断点位置输出加括号的矩阵链。

void Traceback(int i,int j,int s[][N]){//根据输入的i、j限定矩阵链的始末位置，s存储断链点
    if(i==j)       //当单个矩阵时，输出矩阵（递归出口）
    {
        cout<<"A"<<i;
    }
    else       //即i<j 按照最佳断点一分为二，接着继续递归
    {
        cout<<"("; //最外侧的括号
        Traceback(i,s[i][j],s);
        Traceback(s[i][j]+1,j,s);
        cout<<")";
    }
}

完整代码如下：

#include <iostream>
using namespace std;
const int N = 100;
int p[N];//行列值
int m[N][N];//最优解
int s[N][N];//断点位置

void MatrixChain(int n){//构造最优解，p存储数据的维数
    int r, i, j, k;
	for(i=0;i<=n;i++)//初始化对角线
		m[i][i]=0;
	for(r=2;r<=n;r++){//r个矩阵连乘
		for (i=1;i<=n-r+1;i++){
			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;
			for(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;//记录断开位置
				}
			}
		}
	}
}

void Traceback(int i,int j){
    if(i==j){
        cout<<"A"<<i;
    }
    else{
        cout<<"(";
        Traceback(i,s[i][j]);
        Traceback(s[i][j]+1,j);
        cout<<")";
    }
}

int main()
{
    int n;//n个矩阵
	cin >> n;
	int i, j;
	for (i=0;i<=n;i++)
        cin>>p[i];
    //int p[6]={30,35,15,5,10,20,25};
    //int q[5]={5,10,4,6,10,2};
	MatrixChain(n);
	cout<<"ans:";
	Traceback(1,n);
	cout<<endl;
    return 0;
}

终于弄懂了过程！！！真心有点难理解