《算法导论》读书笔记之动态规划—矩阵链乘法

谁借莪１个温暖的怀抱￠ 2021-12-10 16:47 351阅读 0赞

**前言：**今天接着学习动态规划算法，学习如何用动态规划来分析解决矩阵链乘问题。首先回顾一下矩阵乘法运算法，并给出C++语言实现过程。然后采用动态规划算法分析矩阵链乘问题并给出C语言实现过程。

**1、矩阵乘法**

![10103947-1c71326c75c64c26ae52f891fbf8c505.png][]

从定义可以看出： **只有当矩阵*A*的列数与矩阵*B*的行数相等时*A*×*B*才有意义。一个*m*×*r*的矩阵A左乘一个*r*×n的矩阵*B*，会得到一个*m*×*n*的矩阵C**。在计算机中，一个矩阵说穿了就是一个二维数组。一个m行r列的矩阵可以乘以一个r行n列的矩阵，得到的结果是一个m行n列的矩阵，其中的第i行第j列位置上的数等于前一个矩阵第i行上的r个数与后一个矩阵第j列上的r个数对应相乘后所有r个乘积的和。采用C++语言实现完整的两个矩阵乘法，程序如下所示：

![复制代码][copycode.gif]

1 #include <iostream>
     2 using namespace std;
     3 #define A_ROWS        3
     4 #define A_COLUMNS     2
     5 #define B_ROWS        2
     6 #define B_COLUMNS     3
     7 void matrix_multiply(int A[A_ROWS][A_COLUMNS],int B[B_ROWS][B_COLUMNS],int C[A_ROWS][B_COLUMNS]);
     8 int main()
     9 {
    10     int A[A_ROWS][A_COLUMNS] = {1,0,
    11                                 1,2,
    12                                 1,1};
    13     int B[B_ROWS][B_COLUMNS] = {1,1,2,
    14                                 2,1,2};
    15     int C[A_ROWS][B_COLUMNS] = {0};
    16     matrix_multiply(A,B,C);
    17     for(int i=0;i<A_ROWS;i++)
    18     {
    19         for(int j=0;j<B_COLUMNS;j++)
    20             cout<<C[i][j]<<" ";
    21         cout<<endl;
    22     }
    23     return 0;
    24 }
    25 void matrix_multiply(int A[A_ROWS][A_COLUMNS],int B[B_ROWS][B_COLUMNS],int C[A_ROWS][B_COLUMNS])
    26 {
    27     if(A_COLUMNS != B_ROWS)
    28         cout<<"error: incompatible dimensions."<<endl;
    29     else
    30     {
    31         int i,j,k;
    32         for(i=0;i<A_ROWS;i++)
    33             for(j=0;j<B_COLUMNS;j++)
    34             {
    35                 C[i][j] = 0;
    36                 for(k=0;k<A_COLUMNS;k++)
    37                     C[i][j] += A[i][k] * B[k][j]; //将A的每一行的每一列与B的每一列的每一行的乘积求和
    38             }
    39     }
    40 }

![复制代码][copycode.gif]

程序测试结果如下所示：

![10104709-54db72a10dc841ed9f7733753af885dd.png][]

**2、矩阵链乘问题描述**

给定n个矩阵构成的一个链<A1,A2,A3,.......An>，其中i=1,2,...n，矩阵A的维数为pi-1pi，对乘积 A1A2...An 以一种最小化标量乘法次数的方式进行加全部括号。

注意：在矩阵链乘问题中，实际上并没有把矩阵相乘，目的是确定一个具有最小代价的矩阵相乘顺序。找出这样一个结合顺序使得相乘的代价最低。

**3、动态规划分析过程**

**1）最优加全部括号的结构**

动态规划第一步是寻找一个最优的子结构。假设现在要计算AiAi+1....Aj的值，计算Ai...j过程当中肯定会存在某个k值（i<=k<j）将Ai...j分成两部分，使得Ai...j的计算量最小。分成两个子问题Ai...k和Ak+1...j，需要继续递归寻找这两个子问题的最优解。

有分析可以到最优子结构为：**假设AiAi+1....Aj的一个最优加全括号把乘积在Ak和Ak+1之间分开，则Ai..k和Ak+1..j也都是最优加全括号的。**

**2）一个递归解**

设m\[i,j\]为计算机矩阵Ai...j所需的标量乘法运算次数的最小值，对此计算A1..n的最小代价就是m\[1,n\]。现在需要来递归定义m\[i,j\]，分两种情况进行讨论如下：

**当i==j时：m\[i,j\] = 0，（此时只包含一个矩阵）**

**当i<j 时：从步骤1中需要寻找一个k（i≤k＜j）值，使得m\[i,j\] ＝min\{m\[i,k\]+m\[k+1,j\]+pi-1pkpj\} （i≤k＜j）。**

3）计算最优代价

虽然给出了递归解的过程，但是在实现的时候不采用递归实现，而是借助辅助空间，使用自底向上的表格进行实现。设矩阵Ai的维数为pi-1pi，i=1,2.....n。输入序列为：p=<p0，p1,...pn>，length\[p\] = n+1。使用m\[n\]\[n\]保存m\[i,j\]的代价，s\[n\]\[n\]保存计算m\[i,j\]时取得最优代价处k的值，最后可以用s中的记录构造一个最优解。书中给出了计算过程的伪代码，摘录如下：

![复制代码][copycode.gif]

1 MAXTRIX_CHAIN_ORDER(p)
     2   n = length[p]-1;
     3   for i=1 to n
     4       do m[i][i] = 0;
     5   for t = 2 to n  //t is the chain length
     6        do for i=1 to n-t+1
     7                      j=i+t-1;
     8                      m[i][j] = MAXLIMIT;
     9                      for k=i to j-1
    10                             q = m[i][k] + m[k+1][i] + qi-1qkqj;
    11                             if q < m[i][j]
    12                                then m[i][j] = q;
    13                                     s[i][j] = k;
    14   return m and s;

![复制代码][copycode.gif]

**MATRIX\_CHAIN\_ORDER具有循环嵌套，深度为3层，运行时间为O(n3)。如果采用递归进行实现，则需要指数级时间Ω(2n)，因为中间有些重复计算。**递归是完全按照第二步得到的递归公式进行计算，递归实现如下所示：

![复制代码][copycode.gif]

1 int recursive_matrix_chain(int *p,int i,int j,int m[N+1][N+1],int s[N+1][N+1])
     2 {
     3     if(i==j)
     4        m[i][j] = 0;
     5     else
     6     {
     7         int k;
     8         m[i][j] = MAXVALUE;
     9         for(k=i;k<j;k++)
    10         {
    11             int temp = recursive_matrix_chain(p,i,k,m,s) +recursive_matrix_chain(p,k+1,j,m,s) + p[i-1]*p[k]*p[j];
    12             if(temp < m[i][j])
    13             {
    14                 m[i][j] = temp;
    15                 s[i][j] = k;
    16             }
    17         }
    18     }
    19     return m[i][j];
    20 }

![复制代码][copycode.gif]

对递归算计的改进，可以引入备忘录，采用自顶向下的策略，维护一个记录了子问题的表，控制结构像递归算法。完整程序如下所示：

![复制代码][copycode.gif]

1 int memoized_matrix_chain(int *p,int m[N+1][N+1],int s[N+1][N+1])
     2 {
     3     int i,j;
     4     for(i=1;i<=N;++i)
     5         for(j=1;j<=N;++j)
     6         {
     7            m[i][j] = MAXVALUE;
     8         }
     9     return lookup_chain(p,1,N,m,s);
    10 }
    11 
    12 int lookup_chain(int *p,int i,int j,int m[N+1][N+1],int s[N+1][N+1])
    13 {
    14     if(m[i][j] < MAXVALUE)
    15         return m[i][j]; //直接返回，相当于查表
    16     if(i == j)
    17         m[i][j] = 0;
    18     else
    19     {
    20         int k;
    21         for(k=i;k<j;++k)
    22         {
    23             int temp = lookup_chain(p,i,k,m,s)+lookup_chain(p,k+1,j,m,s) + p[i-1]*p[k]*p[j];  //通过递归的形式计算，只计算一次，第二次查表得到
    24             if(temp < m[i][j])
    25             {
    26                 m[i][j] = temp;
    27                 s[i][j] = k;
    28             }
    29         }
    30     }
    31     return m[i][j];
    32 }

![复制代码][copycode.gif]

**4)构造一个最优解**

第三步中已经计算出来最小代价，并保存了相关的记录信息。因此只需对s表格进行递归调用展开既可以得到一个最优解。书中给出了伪代码，摘录如下：

![复制代码][copycode.gif]

1 PRINT_OPTIMAL_PARENS(s,i,j)
    2   if i== j 
    3      then print "Ai"
    4   else
    5      print "(";
    6      PRINT_OPTIMAL_PARENS(s,i,s[i][j]);
    7      PRINT_OPTIMAL_PARENS(s,s[i][j]+1,j);
    8      print")";

![复制代码][copycode.gif]

**4、编程实现**

采用C++语言实现这个过程，现有矩阵A1(30×35)、A2(35×15)、A3(15×5)、A4(5×10)、A5(10×20)、A6(20×25)，得到p=<30,35,15,5,10,20,25>。实现过程定义两个二维数组m和s，为了方便计算其第一行和第一列都忽略，行标和列标都是1开始。完整的程序如下所示：

按 Ctrl+C 复制代码

程序测试结果如下所示：

![10113103-6237493a4ba04322accf6fa507db197c.png][]

**5、总结**

动态规划解决问题关键是分析过程，难度在于如何发现其子问题的结构及子问题的递归解。这个需要多多思考，不是短时间内能明白。在实现过程中遇到问题就是数组，数组的下标问题是个比较麻烦的事情，如何能够过合理的去处理，需要一定的技巧。

From: http://www.cnblogs.com/Anker/archive/2013/03/10/2952475.html

转载于:https://www.cnblogs.com/wushuaiyi/p/4509238.html

[10103947-1c71326c75c64c26ae52f891fbf8c505.png]: /images/20211210/bd902628aa114d138a60928324de6c19.png
[copycode.gif]: /images/20211210/09d3763e1578496aa485dee25a593abb.png
[10104709-54db72a10dc841ed9f7733753af885dd.png]: /images/20211210/ff320ef3a73d4ba69b84b42074d2d049.png
[10113103-6237493a4ba04322accf6fa507db197c.png]: /images/20211210/6554e3756d884346989d0b000c044444.png