【动态规划算法】【Python实现】矩阵连乘问题

文章目录

• • @[toc]
  • • 问题描述
    • 穷举搜索法
    • 动态规划算法
    • • 分析最优解的结构
      • 建立递归关系
      • • $m[i][j]$递归方程
      • 计算最优值
      • 构造最优解
      • `Python`实现
      • 时间复杂性

问题描述

• 给定$n$个矩阵$A_1$，$A_2$，$\cdots$，$A_n$，其中$A_i$与$A_{i+1}$是可乘的（$i=1$，$2$，$\cdots$，$n-1$），$A_i$的维数为$p_{i-1}\times p_i,i=1,2,\cdots ,n$，确定这$n$个矩阵的连乘积$A_1{A}_2\cdots A_n$的计算次序，使得以此计算矩阵连乘积需要的数乘次数最少

穷举搜索法

• 对于$n$个矩阵的连乘积，设有不同的计算次序$P(n)$
• 可以先在第$k$个（$k=1$，$2$，$\cdots$，$n-1$）和第$k+1$个矩阵之间将原矩阵序列分为两个矩阵子序列，然后分别对这两个矩阵子序列完全加括号，最后对所得的结果加括号，得到原矩阵序列的一种完全加括号方式
• $P(n)$递归方程

$$P(n)=\{\begin{array}{ll}1,& n=1\\ \sum _{k=1}^{n-1}P(k)P(n-k),& n>1\end{array}$$

• • $P(n)$实际上是$Catalan$数，即$P(n)=C(n-1)$，$C(n)=\frac{1}{n+1}\left(\begin{array}{c}2n\\ n\end{array}\right)=\Omega (4^n/n^{3/2})$，$P(n)$是随$n$的增长呈指数增长的

动态规划算法

分析最优解的结构

• 考察计算$A[1:n]$的最优计算次序
• 设这个计算次序在矩阵$A_k(1\le k<n)$和$A_{k+1}$之间将矩阵链断开
• 计算$A[1:n]$的最优次序所包含的计算矩阵子链$A[1:k]$和$A[k+1:n]$的次序也是最优的
  • 如果计算$A[1:k]$的次序需要的计算量更少，则用此次序替换原来计算$A[1:k]$的次序，得到的计算$A[1:n]$的计算量将比最优次序所需计算量更少，这是一个矛盾
• 矩阵连乘积计算次序问题的最优解包含其子问题的最优解，这种性质称为最优子结构性质

建立递归关系

• 设计算$A[i:j](1\le i\le j\le n)$所需的最少数乘次数为$m[i][j]$
• 当$i=j$时，$m[i][i]=0(i=1,2,\cdots ,n)$
• 当$i<j$时，若计算$A[i:j]$的最优次序在$A_k(i\le k<j)$和$A_{k+1}$之间断开，则$m[i][j]=m[i][k]+m[k+1][j]+p_{i-1}{p}_k{p}_j$
• 将对应$m[i][j]$的断开位置$k$记为$s[i][j]$

$m[i][j]$递归方程

m [ i ] [ j ] = { 0 , i = j min ⁡ i ≤ k < j {   m [ i ] [ k ] + m [ k + 1 ] [ j ] + p i − 1 p k p j   } , i < j m[i][j] = \begin{cases} 0 , & i = j \\ \min\limits_{i \leq k < j}{\set{m[i][k] + m[k + 1][j] + p_{i - 1} p_{k} p_{j}}} , & i < j \end{cases} m[i][j]=⎩ ⎨ ⎧​0,i≤k<jmin​{m[i][k]+m[k+1][j]+pi−1​pk​pj​},​i=ji<j​

计算最优值

• 对于$1\le i\le j\le n$不同的有序对$(i,j)$对应不同的子问题
• 不同子问题的个数最多只有$\left(\begin{array}{c}n\\ 2\end{array}\right)+n=\theta (n^2)$个

构造最优解

• $s[i][j]$中的数表明，计算矩阵链$A[i:j]$的最佳方式是在矩阵$A_k$和$A_{k+1}$之间断开，即最优加括号方式为$(A[i:k])(A[k+1:j])$

Python实现

def matrix_chain_order(p):
    n = len(p) - 1  # 矩阵的数量

    m = [[float('inf')] * (n + 1) for _ in range(n + 1)]  # 存储最小乘法次数的矩阵
    s = [[0] * (n + 1) for _ in range(n + 1)]  # 存储最佳分割位置的矩阵

    # 对角线上的子问题, 乘法次数为 0
    for i in range(1, n + 1):
        m[i][i] = 0

    # l 表示子问题的规模
    for l in range(2, n + 1):
        for i in range(1, n - l + 2):
            j = i + l - 1

            for k in range(i, j):
                # 计算 m[i][j] 的值
                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, s


def print_optimal_parens(s, i, j):
    if i == j:
        print(f'A{i}', end='')
    else:
        print('(', end='')

        print_optimal_parens(s, i, s[i][j])
        print_optimal_parens(s, s[i][j] + 1, j)

        print(')', end='')


matrix_sizes = [30, 35, 15, 5, 10, 20, 25]

m, s = matrix_chain_order(matrix_sizes)

n = len(matrix_sizes) - 1

print('最佳乘法顺序:', end='')
print_optimal_parens(s, 1, n)
print()

minimum_multiplications = m[1][n]
print('最小乘法次数:', minimum_multiplications)

最佳乘法顺序:((A1(A2A3))((A4A5)A6))
最小乘法次数: 15125

时间复杂性

• 三重循环的总次数为$O(n^3)$，计算时间上界为$O(n^3)$