动态规划之矩阵链乘法

虽然我们在《动态规划之钢条切割》:http://blog.csdn.net/ii1245712564/article/details/44464689里面已经介绍过动态规划的相关知识，这里还是提一下动态的两大特征：

1.最优子结构，原问题的解是在子问题的最优解的基础上建立起来的

2.重叠子问题，即在不断求解子问题的过程中会产生相同的子问题，而不是一直产生新的子问题

1.问题描述

现在我们来看一个例子，假设有三个矩阵：

A1<10,100> A2<100,5> A3<5,50>

令B=A1A2A3,那么求得矩阵B共需要花费多少代价呢？

先看一下这种情况：

A1乘A2的代价为10*100*5=5000,然后得到一个A4<10,5>的矩阵

A4再乘上A3的代价为10*5*50=2500

总共花费代价为5000+2500=7500

我们知道矩阵相乘是不满足乘法交换律，而是满足乘法结合律的，即我们不可以改变矩阵链中各个矩阵的位置，但是我们可以改变矩阵在矩阵链中相乘的先后顺序

那么就有这种情况：

先A2乘A3，代价为1000*5*50=25000,得到矩阵A4<100,50>

A1再乘上A4，代价为10*100*50=50000

总共花费代价为25000+50000=75000

比较上面两种情况，我们知道，通过矩阵的结合律改变矩阵相乘的先后顺序，矩阵相乘的总代价是由区别的。

问题：现在给定一个含有n个矩阵的矩阵链<A1,A2,A3,...,An>（矩阵链输入合法），通过结合律改变矩阵相乘的顺序（通过加括号方法），使得矩阵链花费的代价最小。

2.问题分析

看到这道题，我们首先想到的是用蛮力法来枚举所有的情况，然后找到最小的花费，但是这样做真的可行吗？

我们假设p(n)为n个矩阵链括号化方案的数量：

通过计算，我们知道，通过蛮力法来求解这个问题，那么括号化方案的数量将是与n呈指数增长的关系的。

那我们通过运用动态规划的策略来解决这个问题：

1.刻画一个最优解的结构特征。

2.递归的定义最优解的值

3.计算最优解时，一般采用自下而上的方法来求解

4.利用计算出的信息构造出一个最优解

步骤一：最优括号化方案的结构特征

动态规划的第一步是寻找最优子结构，通过子问题的最优解来构造原问题的最优解。为方便起见，我们用A<i...j>表示矩阵链<Ai,Ai+1,...,Aj-i,Aj>,为了将矩阵链括号化，我们需要在位置k(k>=i && k < j)处将矩阵链分开为A<i...k>和A<k+1...j>。也就是说我们首先计算A<i...k>和A<k+1...j>的子矩阵的代价，然后再加上A<i...k>,A<k+1...j>两者相乘的代价。

步骤二：一个递归解决方案

我们假设：

m[i...j]为矩阵链A<i...j>的标量乘法次数的最小值，那么矩阵链A<1...j>的标量乘法次数为m[1...j]

ri为矩阵Ai的行数，ci为矩阵Ai的列数，那么A1<50,500> r1 = 50,c1=500

于是:

当A<i..j>=Ai,即矩阵链里矩阵数量为1时，m[i...j]=0

当j>i时,首先计算两个子问题的最小乘法次数，分别是m[i...k]和m[k+1...j],然后加上A<i...k>乘上A<k+1...j>的乘法次数ri*ck*cj，有

m[i...j]=m[i...k]+m[k+1...j]+ri*ck*cj

上面的公式中k是已知的，但是在实际中我们是不知道的，于是K有[i...j)种选择，那么递归方程式有：

• m[i...j]=0，i=j
• m[i...j]=min{ m[i,i]+m[i+1,j]+ri*ci*cj , m[i,i+1]+m[i+2,j]+ri*c(i+1)*cj ,  ... , m[i,j-1]+m[j,j]+ri*c(j-1)*cj}

m[i,j]保存的是最小的相乘次数，我们用s[i,j]来保存在A<i...j>里面的括号的分割点

步骤三：计算最优代价

我们知道在矩阵A<1..n>中，矩阵的子问题数量为：，因为每个子问题都是在在1-n之间选择两个不相等的点i,j，即为，最后再加上i,j相等的情况n。而每一个子问题都是从A<i...j>里面选出一个k出来，有n种选择，所以，问题的规模=子问题的数量*子问题有多少种选择=O(n^2)*O(n)=O(n^3)

这里首先采用自上而下带备份的迭代方法实现：

#include <utility>
#include <iostream>
#include <vector>
#include <limits.h>

using namespace std;

#define MAXTRIX_CHAIN_LENGTH 3

unsigned int mulCost[MAXTRIX_CHAIN_LENGTH][MAXTRIX_CHAIN_LENGTH];
unsigned int dealMinCost(std::vector<pair<unsigned ,unsigned> > & matrixChain ,\
                         unsigned left ,\
                         unsigned right);
/**
 * 求解最小花费
 * @param  matrixChain 矩阵链
 * @return             最小花费
 */
unsigned int getMinCost(std::vector<pair<unsigned , unsigned> > & matrixChain)
{
    if(matrixChain.size() <=1)
        return 0;
    dealMinCost(matrixChain ,0, matrixChain.size()-1);
    return mulCost[0][MAXTRIX_CHAIN_LENGTH-1];
}

unsigned int dealMinCost(std::vector<pair<unsigned ,unsigned> > & matrixChain ,\
                         unsigned left ,\
                         unsigned right)
{
    if(left >= right)
        return 0;
    if(mulCost[left][right]!=0)//如果之前算过，那么就直接返回之前的结果
        return mulCost[left][right];
    unsigned int min = UINT_MAX;
    for (int i = left; i < right; ++i)
    {
        unsigned int temp = dealMinCost(matrixChain ,left,i)+\
                            dealMinCost(matrixChain , i+1, right)+\
                            matrixChain[left].first*matrixChain[i].second*matrixChain[right].second;
        if(min > temp)
            min = temp;
    }
    mulCost[left][right] = min;
    return min;
}


int main(int argc, char const *argv[])
{
    vector<pair<unsigned int , unsigned int> > matrixChain;
    matrixChain.push_back(pair<unsigned int ,unsigned int>(10,100));
    matrixChain.push_back(pair<unsigned int ,unsigned int>(100,5));
    matrixChain.push_back(pair<unsigned int ,unsigned int>(5,50));

    cout<<"the min cost is:"<<getMinCost(matrixChain)<<endl;
    return 0;
}

下面这个方法是自下而上的算法：

#include <utility>
#include <iostream>
#include <vector>
#include <limits.h>

using namespace std;

#define MAXTRIX_CHAIN_LENGTH 3

unsigned int mulCost[MAXTRIX_CHAIN_LENGTH][MAXTRIX_CHAIN_LENGTH];
unsigned int dealMinCost(std::vector<pair<unsigned ,unsigned> > & matrixChain ,\
                         unsigned left ,\
                         unsigned right);
/**
 * 求解最小花费
 * @param  matrixChain 矩阵链
 * @return             最小花费
 */
unsigned int getMinCost(std::vector<pair<unsigned , unsigned> > & matrixChain)
{
    if(matrixChain.size() <=1)
        return 0;
    dealMinCost(matrixChain ,0, matrixChain.size()-1);
    return mulCost[0][MAXTRIX_CHAIN_LENGTH-1];
}


unsigned int dealMinCost(std::vector<pair<unsigned ,unsigned> > & matrixChain ,\
                         unsigned left ,\
                         unsigned right)
{
    for (int n = 0; n < right-left+1  ; ++n)
    {
        int i,j;
        for (i = left , j = left+n ; i <= right-n && j <= right ; ++i , ++j)
        {
            if(i == j)
            {
                mulCost[i][j] = 0;
                continue;
            }
            unsigned int min = UINT_MAX;
            for (int k = i; k < j; ++k)
            {
                unsigned int temp = mulCost[i][k]+mulCost[k+1][j]+\
                                    matrixChain[i].first*matrixChain[k].second*matrixChain[j].second;
                if(min > temp)
                    min = temp;
            }
            mulCost[i][j] = min;
        }
    }
}


int main(int argc, char const *argv[])
{
    vector<pair<unsigned int , unsigned int> > matrixChain;
    matrixChain.push_back(pair<unsigned int ,unsigned int>(10,100));
    matrixChain.push_back(pair<unsigned int ,unsigned int>(100,5));
    matrixChain.push_back(pair<unsigned int ,unsigned int>(5,50));

    cout<<"the min cost is:"<<getMinCost(matrixChain)<<endl;
    return 0;
}

步骤四：构造最优解

在A<i...j>中，我们需要对A<i...j>进行一次划分，以得到最优解，所以我们只需记录A<i...j>的划分点k即可，所以代码如下：

其中solution[i][j]就是用来保存A<i...j>的划分点k的，printSolution是用来递归打印解决方案的函数

#include <utility>
#include <iostream>
#include <vector>
#include <limits.h>

using namespace std;

#define MAXTRIX_CHAIN_LENGTH 6

unsigned int mulCost[MAXTRIX_CHAIN_LENGTH][MAXTRIX_CHAIN_LENGTH];
unsigned int solution[MAXTRIX_CHAIN_LENGTH][MAXTRIX_CHAIN_LENGTH];
unsigned int dealMinCost(std::vector<pair<unsigned ,unsigned> > & matrixChain ,\
                         unsigned left ,\
                         unsigned right);
/**
 * 求解最小花费
 * @param  matrixChain 矩阵链
 * @return             最小花费
 */
unsigned int getMinCost(std::vector<pair<unsigned , unsigned> > & matrixChain)
{
    if(matrixChain.size() <=1)
        return 0;
    dealMinCost(matrixChain ,0, matrixChain.size()-1);
    return mulCost[0][MAXTRIX_CHAIN_LENGTH-1];
}


unsigned int dealMinCost(std::vector<pair<unsigned ,unsigned> > & matrixChain ,\
                         unsigned left ,\
                         unsigned right)
{
    for (int n = 0; n < right-left+1  ; ++n)
    {
        int i,j;
        for (i = left , j = left+n ; i <= right-n && j <= right ; ++i , ++j)
        {
            if(i == j)
            {
                mulCost[i][j] = 0;
                continue;
            }
            unsigned int min = UINT_MAX;
            unsigned int pos;
            for (int k = i; k < j; ++k)
            {
                unsigned int temp = mulCost[i][k]+mulCost[k+1][j]+\
                                    matrixChain[i].first*matrixChain[k].second*matrixChain[j].second;
                if(min > temp)
                {
                    min = temp;
                    pos = k;
                }
            }
            mulCost[i][j] = min;
            solution[i][j] = pos;
        }
    }
}


void printSolution(int left , int right)
{
    if(left+2 > right)
        return;
    cout<<"from "<<left<<" to "<<right<<" the bracket pos is "<<solution[left][right]<<endl;
    printSolution(left,solution[left][right]);
    printSolution(solution[left][right]+1,right);
}

int main(int argc, char const *argv[])
{
    vector<pair<unsigned int , unsigned int> > matrixChain;
    matrixChain.push_back(pair<unsigned int ,unsigned int>(30,35));
    matrixChain.push_back(pair<unsigned int ,unsigned int>(35,15));
    matrixChain.push_back(pair<unsigned int ,unsigned int>(15,5));
    matrixChain.push_back(pair<unsigned int ,unsigned int>(5,10));
    matrixChain.push_back(pair<unsigned int ,unsigned int>(10,20));
    matrixChain.push_back(pair<unsigned int ,unsigned int>(20,25));

    cout<<"the min cost is:"<<getMinCost(matrixChain)<<endl;
    printSolution(0,MAXTRIX_CHAIN_LENGTH-1);
    while(1);
    return 0;
}