矩阵链乘——动态规划

普通矩阵相乘算法

//普通矩阵相乘
#define ROW_A 2
#define COL_A 2
#define COL_B 3
void matrix_mul(int matA[ROW_A][COL_A], int matB[2][COL_B], int C[ROW_A][COL_B])
{

    for (int i = 0; i < ROW_A; i++) {
        for (int j = 0; j < COL_B; j++) {
            //A的行必须等于B的列
            for (int k = 0; k < COL_A; k++) {
                C[i][j] += (matA[i][k] * matB[k][j]);
            }
        }
    }
}

两个相容矩阵相乘，相容是指矩阵A的行必须等于矩阵B的列。

利用矩阵相乘的定义计算，Cij = ∑ Aik * Bkj （k = 1……A.col）。

需要注意的是Aik 和 Bhj 规模的两个矩阵相乘后，矩阵规模会变为 Cij，而两个矩阵相乘的代价是i * j* k.

之后的所有算法都会用到矩阵的规模和代价。

递归求矩阵链乘算法

约定 Ai……j 表示矩阵Ai 一直乘到 Aj 的一个矩阵链

首先发掘最优子结构。

1. 做出一个选择，选择n个矩阵在第k个位置划分
2. 假定我已经知道在第k个位置划分后代价最优
3. 给定可获得最优解的选择后，确定选择会产生那些子问题，（产生Ai……k 和 Ak+1……..j 两个矩阵链乘的子问题），以及如何刻画子问题解空间
4. 剪切粘贴技术证明，原问题的最优解包含子问题的最优解。因为如果子问题的解不是最优解，那么把子问题的最优解代入原问题中就能得到更优的原问题解，所以证明了子问题的解也必定是最优的。
5. 原问题产生的两个子问题解空间形式不太一样，例如A1……k ，Ak+1…….j ，后者的形式和原问题A1…….j的形式不同，必须允许子问题i，j均可变，所以用m[i][j]记录解空间.

最优子结构的体现

原问题的最优解涉及多少个子问题？

在确定最优解使用那些子问题时，需要考察多少种选择？

答案是：最优解使用两个子问题，需要考察j-i种选择。

两个子问题是Ai……k ,和Ak+1…….j 各自的最优解，还要考察k在i 到 j 之间取值最优的。

用子问题的解递归定义原问题的最优解

m[i][j] 表示Ai….j个矩阵相乘的最小代价，最优解就是m[1][n]——即从第一个矩阵乘到最后一个矩阵的代价。

• 递归的定义m[i][j]

i = j 时，就是一个矩阵乘自己，代价为0

i < j 时，在k处划分并且考察所有k的选择，m[i][j] = min（m[i][k] +m[k+1][j] + pi-1 * pk *pj）(k = i…….j-1)

• 在本题中因为给出的矩阵代价是p数组且有规律，所以可以用p表达出来前后两个矩阵相乘代价

上面是求解动归的最重要部分，基本可以直接写出递归的方法了↓

//矩阵规模，例如A1=30*35，A2=35*15
int p[] = { 30,35,15,5,10,20,25 };
#define N  6        //矩阵长度
int m[N+1][N+1] = { 0 };
int s[N + 1][N + 1] = {0};
//直接递归法,i=1,j=6
int recursive(int p[], int i, int j)
{
    if (i == j) {
        return 0;
    }
    m[i][j] = _CRT_INT_MAX;
    for (int k = i; k < j; k++) {
        int q = recursive(p, i, k) + recursive(p, k + 1, j) + p[i - 1] * p[k] * p[j];
        if (q < m[i][j]) {
            m[i][j] = q;
        }
    }
    return m[i][j];
}

带备忘的递归算法

当递归调用第一次遇到子问题就求其解，记录。以后再遇到直接查表返回。

//带备忘的递归
void memorize(int m[N + 1][N + 1]) {
    for (int i = 0; i < N + 1; i++) {
        for (int j = 0; j < N + 1; j++) {
            m[i][j] = _CRT_INT_MAX;
        }
    }
}
int lookup(int p[], int i, int j)
{
    //查表
    if (m[i][j] < _CRT_INT_MAX) {
        return m[i][j];
    }
    if (i == j) {
        return 0;
    }
    for (int k = i; k < j; k++) {
        int q = recursive(p, i, k) + recursive(p, k + 1, j) + p[i - 1] * p[k] * p[j];
        if (q < m[i][j]) {
            m[i][j] = q;
        }
    }

    return m[i][j];
}

自底向上

为了实现自底向上法，必须确定计算m[i][j]时访问哪些表项。由上面的递归定义式看，m[i][j] 的一次计算需要求所有k的选择的子问题的解中最优的，而k在i…j中取值，意味着需要访问的子问题的矩阵规模一定小于原问题。所以设计求解顺序为按矩阵链长度从短到长。

而且自底向上需要求出所有的子问题解，m[i][j]合法解空间会被填满，所以循环考虑每个长度，每个起始结束点的矩阵链，最后返回最优代价是m[1][N].

//自底向上
int matrix_order(int p[])
{
    //数组下标一律和矩阵序号对齐
    for (int k = 1; k <= N; k++) {
        m[k][k] = 0;
    }
    //设计求解顺序，按长度
    for (int len = 2; len <= N; len++) {
        //i开始矩阵下标，j和i间隔len个矩阵
        for (int i = 1; i <= N - len+1; i++) {
            int j = i + len - 1;
            m[i][j] = _CRT_INT_MAX;
            for (int k = i; k < j; k++) {
                int q = m[i][k] + m[k+1][j] + p[i - 1] * p[k] * p[j];
                if (q < m[i][j]) {
                    m[i][j] = q;
                    s[i][j] = k;
                }
            }
        }

    }
    return m[1][N];
}
//构造最优解
void print(int s[N+1][N+1], int i, int j)
{
    if (i == j) {
        printf("A%d", i);
    }
    else {
        printf("(");
        print(s, i, s[i][j]);
        print(s, s[i][j]+1, j);
        printf(")");
    }
}