区间DP | 1:矩阵链乘问题(含优化) —— 例题:矩阵链乘、合并石子

矩阵链乘法问题: 给定 n 个矩阵的链 <A1, A2, A3, ..., An>,矩阵 Ai 的规模为 p_{i-1} \cdot p_{i} (0\leq i \leq n) 。求完全括号化方案,使得计算乘积A1A2...An所需要标量的乘法次数最少。


目录

一、算法剖析

1、最优括号化方案的结构特征

2、 一个递归的求解方案

二、算法实现

1、计算最优代价

暂未优化的方法 —— 时间复杂度为 O(n^3)

优化后的方法 —— 时间复杂度降至 O(n^2)

2、构造最优解

三、例题

1、矩阵链乘问题

完整AC代码1:(未优化的时间复杂度 O(n^3) 

完整AC代码2:(优化的时间复杂度 O(n^2) ) 

2、合并石子

完整AC代码1:(未优化的时间复杂度O(n^3) 

完整AC代码2:(优化的时间复杂度 O(n^2) ) 

end 



两个矩阵 A 和 B 只有相容,即 A 的列数等于 B 的行数时,才能相乘。如果 A 是 p×q 矩阵 ,B 是 q×r 矩阵,那么乘积 C 是 p× r 矩阵。计算 C 的进行了 pqr 次乘法运算。我们采用乘法进行的次数来表示计算的代价。 

对于多个矩阵连乘的问题,我们称作矩阵链乘问题。矩阵的乘法满足结合律,因此任何加括号的方法并不会影响多个矩阵相乘的结果。不同的加括号的方案(改变计算顺序)会影响矩阵链乘法计算的代价:


一、算法剖析

1、最优括号化方案的结构特征

动态规划方法的第一步是寻找最优子结构,然后就可以利用这种子结构,从子问题的最优解构造出原问题的最优解。

(为方便起见,我们用 A_{i..j } 表示A_iA_{i+1}...A_j的乘积结果矩阵)

  • 整体往局部看:对于每一个矩阵链 A_{i..j},如果对其括号化,我们应该在某两个矩阵 A_dA_{d+1} 之间划分开,然后我们再讨论 A_{i..d}A_{d+1..j}的括号化。
  • 再从局部往整体看:计算A_{i..j } 的代价 = 计算 A_{i..d} 的代价 + 计算 A_{d+1..j } 的代价 + 计算 A_{i..d} \cdot A_{d+1..j} 的代价。

因此,为了构造一个矩阵链 A_{i..n } 乘法问题的最优解,我们可以将问题分解成两个子问题: A_{i..d} 和 A_{d+1..n } 最优化括号问题。求出子问题的最优解后,再将子问题的最优解组合起来


2、 一个递归的求解方案

  • 矩阵的规模对应存储在一维数组 p[0..n] 中,矩阵 A_i 的规模为 p_{i-1}\cdot p_{i}
  • 计算的代价对应储存在二维数组 ans[1..n, 1..n] 中,我们用 ans[i, j] 表示矩阵 A_{i..j} 的最小计算代价,那么原问题的最小计算代价最终就在 ans[1, n] 内。 
  • 最优分割点对应储存在二维数组 divide[1..n, 1..n] 中,我们用 divide[i, j] 记录最优代价 ans[i, j] 对应的最优分割点 d。(此数组在打印最优结构时十分有用)

根据上述最优括号化方案的结构特征,我们可以列出状态转移方程:

  • 当 i = j 时: \bg_white ans[i, j] = 0
  • 当 i < j 时: ans[i, j] = min_{i\leq k< j} (ans[i, k] +ans[k+1,j] + p_{i-1}p_kp_j)

ps:在定义时将 ans数组全部初始化为0。


二、算法实现

1、计算最优代价

从算法剖析中已知,本问题需要自底向上的方法:先处理好子矩阵链的最优问题,再总合起来。

那么我们实现时:应该从长度为2的矩阵链开始讨论,然后再在此基础上讨论长度为3的矩阵链...最终讨论长度的n的矩阵链,即得到答案。在某个已确定长度的讨论中,不可以漏去情况,不然会对最优值有影响。


暂未优化的方法 —— 时间复杂度为 O(n^3)

那么我们的循环结构应该是这样:

  • 第一层循环:长度 l = 2..n
  • 第二层循环:讨论每个长度为 l 的矩阵链 A_{i..j}:i = 1..n - l +1, j = i + l + 1
  • 第三层循环:确定最优分割点 k = i..j-1
/* 计算n元矩阵链p的最优代价
 * 动态规划的每一步结果储存在 ans与 divide中 */
void MatrixChainOrder(int p[], int n) {
    /* 矩阵链长度为2到n */
    for (int l = 2; l <= n; l++) {
        /* 讨论长度为l的矩阵链A[i..j] */
        for (int i = 1; i <= n - l + 1; i++) {
            int j = i + l - 1;
            ans[i][j] = INT_MAX;
            /* 依次讨论每一个分割点d:将矩阵链A[i..j]分成A[i..d]和 A[d+1..j] */
            for (int temp, d = i; d < j; d++) {
                temp = ans[i][d] + ans[d + 1][j] + p[i - 1] * p[d] * p[j];
                /* 记录下矩阵链A[i..j]最小的情况 */
                if (temp < ans[i][j]) {
                    ans[i][j] = temp;
                    divide[i][j] = d;
                }
            }
        }
    }
}

优化后的方法 —— 时间复杂度降至O(n^2) 

上面基本方法的三层循环不可避免,但是可大大削减第三层循环的遍历次数,使得时间复杂度降至O(n^2) 。根据矩阵链乘问题模型的特点,其动态规划求解可以使用四边形不等式优化。(链接里有详细的优化原理、时间复杂度的证明过程)

最佳决策点满足如此关系:divide[i][j-1]\leq divide[i][j]\leq divide[i+1][j],所以我们的第三层循环可以优化为:在特定区间 [ divide[ i ][ j -1], divide[i +1][ j ] ]内遍历,其遍历范围远远小于n。

ps:要注意边界情况,需要预先将divide的初值准备好。

/* 计算n元矩阵链p的最优代价
 * 动态规划的每一步结果储存在 ans与 divide中 */
void MatrixChainOrder(int p[], int n) {
    /* 初始化divide数组 */
    for (int i = 1; i <= n; i++)
        divide[i][i] = i;
    /* 矩阵链长度为2到n */
    for (int l = 2; l <= n; l++) {
        /* 讨论长度为l的矩阵链A[i..j] */
        for (int i = 1; i <= n - l + 1; i++) {
            int j = i + l - 1;
            ans[i][j] = INT_MAX;
            /* 在特定区间内:依次讨论每一个分割点d,将矩阵链A[i..j]分成A[i..d]和 A[d+1..j] */
            for (int temp, d = divide[i][j - 1]; d <= divide[i + 1][j]; d++) {
                temp = ans[i][d] + ans[d + 1][j] + p[i - 1] * p[d] * p[j];
                /* 记录下矩阵链A[i..j]最小的情况 */
                if (temp < ans[i][j]) {
                    ans[i][j] = temp;
                    divide[i][j] = d;
                }
            }
        }
    }
}

2、构造最优解

我们在二维数组 divide[i, j] 中已经记录了 A_{i..j} 的最优分界点,由于矩阵链的最优划分可以分解为子链的最优划分,故我们可以递归地打印出来。

/* 根据divide数组,输出A[i..j]的最优情况 */
void PrintDivide(int i, int j) {
    if (i == j)
        printf("A%d", i);
    else {
        int d = divide[i][j];  //找到分界点
        printf("(");
        PrintDivide(i, d);
        PrintDivide(d + 1, j);
        printf(")");
    }
}

三、例题

1、矩阵链乘问题

成绩10开启时间2020年03月10日 星期二 07:55
折扣0.8折扣时间2020年04月7日 星期二 23:55
允许迟交关闭时间2020年04月7日 星期二 23:55

输入:

共两行

第一行 N ( 1<=N<=100 ),代表矩阵个数。

第二行有 N+1 个数,分别为 A1 、 A2 ...... An+1 ( 1<=Ak<=2000 ), Ak 和 Ak+1 代表第 k 个矩阵是个 Ak X Ak+1 形的。

输出:

共两行

第一行 M ,为最优代价。注:测试用例中 M 值保证小于 2^31

第二行为最优顺序。如 (A1((A2A3)A4)) ,最外层也加括号。

注意:测试用例已经保证了输出结果唯一,所以没有AAA的情况.

 测试输入期待的输出时间限制内存限制额外进程
测试用例 1 
  1. 6↵
  2. 30 35 15 5 10 20 25↵
 
  1. 15125↵
  2. ((A1(A2A3))((A4A5)A6))↵
1秒64M

0

 


具体方法上面已经详细讲解啦~

主要注意:当矩阵链中矩阵数目为1的输出要加上括号,不然会wa一个。

完整AC代码1:(未优化的O(n^3)时间复杂度)

//
// Created by LittleCat on 2020/3/10.
//

#include <cstdio>
#include <climits>

using namespace std;
#define MAX 2010

int ans[MAX][MAX] = {0}; // 保存矩阵链 A[i..j]的最小代价
int divide[MAX][MAX] = {0};   // 记录最小代价ans[i,j]对应的分割点

/* 计算n元矩阵链p的最优代价
 * 动态规划的每一步结果储存在 ans与 divide中 */
void MatrixChainOrder(int p[], int n) {
    /* 矩阵链长度为2到n */
    for (int l = 2; l <= n; l++) {
        /* 讨论长度为l的矩阵链A[i..j] */
        for (int i = 1; i <= n - l + 1; i++) {
            int j = i + l - 1;
            ans[i][j] = INT_MAX;
            /* 依次讨论每一个分割点d:将矩阵链A[i..j]分成A[i..d]和 A[d+1..j] */
            for (int temp, d = i; d < j; d++) {
                temp = ans[i][d] + ans[d + 1][j] + p[i - 1] * p[d] * p[j];
                /* 记录下矩阵链A[i..j]最小的情况 */
                if (temp < ans[i][j]) {
                    ans[i][j] = temp;
                    divide[i][j] = d;
                }
            }
        }
    }
}

/* 根据divide数组,输出A[i..j]的最优情况 */
void PrintDivide(int i, int j) {
    if (i == j)
        printf("A%d", i);
    else {
        int d = divide[i][j];  //找到分界点
        printf("(");
        PrintDivide(i, d);
        PrintDivide(d + 1, j);
        printf(")");
    }
}

int main() {
    int n;
    scanf("%d", &n);
    int p[MAX];
    for (int temp, i = 0; i <= n; i++)
        scanf("%d", &p[i]);

    MatrixChainOrder(p, n);
    printf("%d\n", ans[1][n]);
    if (n == 1)
        printf("(A1)");
    else
        PrintDivide(1, n);
    printf("\n");
}


完整AC代码2:(优化的O(n^2)时间复杂度) 

//
// Created by LittleCat on 2020/3/10.
//

#include <cstdio>
#include <climits>

using namespace std;
#define MAX 2010

int ans[MAX][MAX] = {0}; // 保存矩阵链 A[i..j]的最小代价
int divide[MAX][MAX] = {0};   // 记录最小代价ans[i,j]对应的分割点

/* 计算n元矩阵链p的最优代价
 * 动态规划的每一步结果储存在 ans与 divide中 */
void MatrixChainOrder(int p[], int n) {
    /* 初始化divide数组 */
    for (int i = 1; i <= n; i++)
        divide[i][i] = i;
    /* 矩阵链长度为2到n */
    for (int l = 2; l <= n; l++) {
        /* 讨论长度为l的矩阵链A[i..j] */
        for (int i = 1; i <= n - l + 1; i++) {
            int j = i + l - 1;
            ans[i][j] = INT_MAX;
            /* 在特定区间内:依次讨论每一个分割点d,将矩阵链A[i..j]分成A[i..d]和 A[d+1..j] */
            for (int temp, d = divide[i][j - 1]; d <= divide[i + 1][j]; d++) {
                temp = ans[i][d] + ans[d + 1][j] + p[i - 1] * p[d] * p[j];
                /* 记录下矩阵链A[i..j]最小的情况 */
                if (temp < ans[i][j]) {
                    ans[i][j] = temp;
                    divide[i][j] = d;
                }
            }
        }
    }
}

/* 根据divide数组,输出A[i..j]的最优情况 */
void PrintDivide(int i, int j) {
    if (i == j)
        printf("A%d", i);
    else {
        int d = divide[i][j];  //找到分界点
        printf("(");
        PrintDivide(i, d);
        PrintDivide(d + 1, j);
        printf(")");
    }
}

int main() {
    int n;
    scanf("%d", &n);
    int p[MAX];
    for (int temp, i = 0; i <= n; i++)
        scanf("%d", &p[i]);

    MatrixChainOrder(p, n);
    printf("%d\n", ans[1][n]);
    if (n == 1)
        printf("(A1)");
    else
        PrintDivide(1, n);
    printf("\n");
}



2、合并石子

合并石子问题:线性排列着N堆石子,现在要将石子合并成一堆。规定如下:每次只能将相邻的两堆石子合并,合并两堆石子所花费的时间为两堆石子的数量和。求将N堆石子合并成一堆最小花费的时间。(石子分为n堆,石子的数量存储在数组p_{0..n-1}中)

 下面是一个 n = 4 、p = {4, 2, 3, 4} 的例子,最小花费的时间为26。


此题本质就是一个矩阵链乘问题 —— 在链中选择每相邻两节点合并,直到全部合并完毕为止。合并石子的过程也就是在矩阵链中加括号的过程。

此题也是一个具有最优子结构的问题:

为了构造一堆石子 p_{0..n-1} 合并的最优解,我们可以将问题分解成两个子问题: \bg_white p_{i..d} 和 p_{d+1..n }的石子合并问题。求出子问题的最优解后,再将子问题的最优解组合起来

寻找一个递归的动态求解方案:

合并的代价对应储存在二维数组 ans[0..n-1, 0..n-1] 中,我们用 ans[i, j] 表示合并石子堆p_{i..j} 的最小计时间代价,那么原问题的最小计算代价最终就在 ans[1, n] 内。 

状态转移方程:

  • 当 i = j 时: \bg_white ans[i, j] = 0
  • 当 i < j 时: ans[i][j] = min_{i\leq k< j}(ans[i][k] + ans[k+1][j] + \sum_{i\leq t\leq j} p[t])

 ps:在定义时将 ans数组全部初始化为 0。


完整AC代码1:(未优化的O(n^3)时间复杂度)

核心部分是三层循环

  • 第一层循环:长度 l = 2..n
  • 第二层循环:讨论每个长度为 l 的石子堆 p_{i..j}:i = 1..n - l +1, j = i + l + 1
  • 第三层循环:确定最优分割点 k = i..j-1

可见,此代码的时间复杂度为 O(n^3)

#define N 100

#include <cstdio>
#include <climits>

/* 合并石子问题:线性排列着N堆石子,现在要将石子合并成一堆。
 * 规定如下:每次只能将相邻的两堆石子合并,合并两堆石子所花费的时间为两堆石子的数量和。
 * 求将N堆石子合并成一堆最小花费的时间。(石子分为n堆,石子的数量存储在数组p[0..n-1]中)*/
int ans[N][N] = {0};

int MergeStone(int p[], int n) {

    /* 石子堆的个数:从1到n */
    for (int l = 2; l <= n; l++) {
        /* 讨论l个石子的石子堆 p[i..j] */
        for (int i = 0; i < n - l + 1; i++) {
            int j = i + l - 1, sum = 0;
            /* 计算p[i..j]的石子总和 */
            for (int t = i; t <= j; t++)
                sum += p[t];
            ans[i][j] = INT_MAX;
           /* 依次讨论每一个分割点d:将石子堆p[i..j]分成p[i..k]和 A[k+1..j] */
            for (int temp, k = i; k < j; k++) {
                temp = ans[i][k] + ans[k + 1][j] + sum;
                if (temp < ans[i][j])
                    ans[i][j] = temp;
            }
        }
    }
    return ans[0][n - 1];
}

完整AC代码2:(优化的O(n^2)时间复杂度) 

#define N 100

#include <cstdio>
#include <climits>

/* 合并石子问题:线性排列着N堆石子,现在要将石子合并成一堆。
 * 规定如下:每次只能将相邻的两堆石子合并,合并两堆石子所花费的时间为两堆石子的数量和。
 * 求将N堆石子合并成一堆最小花费的时间。(石子分为n堆,石子的数量存储在数组p[0..n-1]中)*/

//优化版本!!!!

int ans[N][N] = {0};
int divide[N][N];

int MergeStone(int p[], int n) {
    /* divide数组初始化*/
    for (int i = 0; i < n; i++)
        divide[i][i] = i;

    /* 石子堆的个数:从1到n */
    for (int l = 2; l <= n; l++) {
        /* 讨论l个石子的石子堆 p[i..j] */
        for (int i = 0; i < n - l + 1; i++) {
            int j = i + l - 1, sum = 0;
            /* 计算p[i..j]的石子总和 */
            for (int t = i; t <= j; t++)
                sum += p[t];
            ans[i][j] = INT_MAX;
            /* 依次讨论特定区间内分割点d:将石子堆p[i..j]分成p[i..k]和 A[k+1..j] */
            for (int temp, k = divide[i][j - 1]; k <= divide[i + 1][j]; k++) {
                temp = ans[i][k] + ans[k + 1][j] + sum;
                if (temp < ans[i][j]) {
                    ans[i][j] = temp;
                    divide[i][j] = k;  //更新p[i..j]的石子堆的最优分割位置
                }

            }
        }
    }
    return ans[0][n - 1];
}

有任何问题欢迎评论交流,如果本文对您有帮助不妨点点赞,嘻嘻~  



end 

欢迎关注个人公众号 鸡翅编程 ”,这里是认真且乖巧的码农一枚。

---- 做最乖巧的博客er,做最扎实的程序员 ----

旨在用心写好每一篇文章,平常会把笔记汇总成推送更新~

在这里插入图片描述

  • 9
    点赞
  • 29
    收藏
    觉得还不错? 一键收藏
  • 2
    评论

“相关推荐”对你有帮助么?

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论 2
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值