第3章 动态规划——矩阵连乘最优计算方式查找

  问题:如何得到n个矩阵连乘A_{p_1p_2}A_{p_2p_3}A_{p_3p_4}...A_{p_{n-1}p_n}A_{p_{n}p_{n+1}}的最少计算次数的计算顺序?先计算A_{p_1p_2}A_{p_2p_3},还是先计算A_{p_2p_3}A_{p_3p_4}?其中,p_k为矩阵的维度。

1、两个矩阵相乘A_{p_1p_2}A_{p_2p_3}

  两个矩阵相乘A_{p_1p_2}A_{p_2p_3}需要多少次运算呢?需要做p_1\times p_2\times p_3次乘法。具体来说,它将得到矩阵(p_1\times p_3)大小的矩阵,该矩阵中每个元素由p_2次乘法得到。

2、多个矩阵连乘

  多个矩阵连乘会有多少种计算方式呢?

  1. 穷举法
    A_{p_1p_2}A_{p_2p_3}A_{p_3p_4}...A_{p_{n-1}p_n}A_{p_{n}p_{n+1}},假设从第k个矩阵断开,分为
    (A_{p_1p_2}A_{p_2p_3}...A_{p_kp_{k+1}})(A_{p_{k+1}p_{k+2}}...A_{p_{n-1}p_n}A_{p_{n}p_{n+1}}),先计算前半部分,再计算后半部分,最后前后两个矩阵相乘。则前半部分有P(k)种计算方式,后半部分P(n-k)种,总共P(k)P(n-k)种计算方式。
     k从1到n-1,则总共有P(n)种:
    P(n)=\sum_{k=1}^{n-1}{P(k)P(n-k)}, P(1)=1
     据说P(n)=C(n-1)=\frac{1}{n} \left( \begin{matrix} 2(n-1) \\ (n-1) \end{matrix} \right) = \Omega \left(4^{n-1}/(n-1)^{3/2}\right)

  2. 上述方法中,会涉及到递归,同时有很多重复运算。
      动态规划就是把重复计算的结果保存下来,下次遇到时只需要查找到已经计算好的结果。即上述穷举过程中,从i到j个矩阵相乘A_{p_{i}p_{i+1}}A_{p_{i+1}p_{i+2}}...A_{p_{j}p_{j+1}},记为A[i:j],会出现很多次,记录下来后,就不用重复计算。去除重复计算后,整个连乘只需要计算n^2次。

  3. 程序
      以A_{30\times 35}A_{35\times 15}A_{15\times 5}A_{5\times 10}A_{10\times 20}A_{20\times 25}矩阵连乘为例,优化后递归实现方式的Python程序如下:

# 输入 一组矩阵的shape,A[[p0, p1], [p1, p2], [p2, p3],  ..., [p_n-1, p_n]]
# 中间存储 A[i,:]到A[j,:]到最小计算次数m[i,j]、最优断开点s
# 输出  A[0,:]到A[n-1,:]到最小计算次数、最优断开点s

import numpy as np 
A = [[30,35], [35,15], [15,5], [5,10], [10,20], [20,25]]
n = np.shape(A)[0] #获得连乘矩阵的个数n
m = np.zeros((n, n, 2), dtype=int) #存储[i,j]从Ai到Aj的[最小计算次数, 最优断开点s]

# A里面相邻矩阵:前一个矩阵列数必须等于后一个矩阵的行数,否则返回错误
def isAillegal(A, i, j):
    # 略
    return False

# 主要函数matrixChain
# 输入  矩阵链A
#      i,j表示从Ai到Aj的 
# 返回  从Ai到Aj的[最小计算次数, 最优断开点s],并记录在m中
def matrixChain(A, i, j):
    # 输入合法性判断
    if i==j:
        return [0,0] #i==j时,m[i,j]的最小计算次数为0
    if i>j or j>=n or i<0:
        return [-1,-1] 
    if isAillegal(A, i, j):
        return [-1,-1]

    # 查找m是否已经有值
    if m[i,j,0]!=0 :
        return m[i,j]

    #相邻矩阵相乘,例如A(P1,P2)*A(P2,P3)需要P1*P2*P3次乘法操作
    if j==(i+1):
        m[i][j][0] = A[i][0]*A[i][1]*A[j][1]
        m[i][j][1] = i
        return m[i,j]

    # 矩阵序列从s处分成前后两部分,循环判断s的位置使得矩阵连乘计算次数最少
    minCount = 0
    minS = i
    for s in range(i,j):
        first = matrixChain(A,i,s)
        last = matrixChain(A,s+1,j)
        count = first[0] + last[0] + A[i][0]*A[s][1]*A[j][1]
        if (s==i) or (count<minCount):
            minCount = count
            minS = s
    m[i][j] = [minCount, minS]
    return m[i,j]

# 输出Ai到Aj的最优计算方式,类似(A0(A1A2))((A3A4)A5))
def traceS(mm, i, j):
    if i==j:
        return 'A'+str(i)  #只剩余一个矩阵,则返回‘Ai’
    if (i+1)==j:
        return 'A'+str(i)+'A'+str(j)+'' #剩余2个矩阵,则返回‘AiAj’
    
    # 否则,分成左右两个部分
    leftStr = '('+traceS(mm, i, mm[i,j,1])+')'
    RightStr = '('+traceS(mm, mm[i,j,1]+1, j)+')'
    return leftStr+RightStr
    
matrixChain(A, 0, n-1)
print(m[:,:,1])
print(traceS(m, 0, n-1))

结果如下:

最优计算次数:
[[    0 15750  7875  9375 11875 15125]
 [    0     0  2625  4375  7125 10500]
 [    0     0     0   750  2500  5375]
 [    0     0     0     0  1000  3500]
 [    0     0     0     0     0  5000]
 [    0     0     0     0     0     0]]
最优断开处:
[[0 0 0 2 2 2]
 [0 0 1 2 2 2]
 [0 0 0 2 2 2]
 [0 0 0 0 3 4]
 [0 0 0 0 0 4]
 [0 0 0 0 0 0]]
最优计算方式:
((A0)(A1A2))((A3A4)(A5))

  直接动态优化算法的Python程序如下:

import numpy as np 
A = [[30,35], [35,15], [15,5], [5,10], [10,20], [20,25]]
n = np.shape(A)[0] #获得连乘矩阵的个数n
m = np.zeros((n, n, 2), dtype=int) #存储[i,j]从Ai到Aj的[最小计算次数, 最优断开点s]

def matrixChainDirect(A):
    for k in range(1, n):  #从2个矩阵连乘开始,计算到n个矩阵连乘。
        for i in range(0,n-k): #从第i个矩阵开始连乘k个矩阵
            minCount = 0
            minS = i
            for s in range(i,i+k): # A_i * A_i+1 *... *A_i+k-1 矩阵连乘从s处断开
                count = m[i][s][0] + m[s+1][i+k][0] + A[i][0]*A[s][1]*A[i+k][1]
                if (s==i) or (count<minCount):
                    minCount = count
                    minS = s
            m[i][i+k] = [minCount, minS]
    return m

matrixChainDirect(A)
print('最优计算次数:')
print(m[:,:,0])
print('最优断开处:')
print(m[:,:,1])
print('最优计算方式:')
print(traceS(m, 0, n-1))

  结果同递归调用方式。

©著作权归作者所有,转载或内容合作请联系作者
  • 序言:七十年代末,一起剥皮案震惊了整个滨河市,随后出现的几起案子,更是在滨河造成了极大的恐慌,老刑警刘岩,带你破解...
    沈念sama阅读 199,711评论 5 468
  • 序言:滨河连续发生了三起死亡事件,死亡现场离奇诡异,居然都是意外死亡,警方通过查阅死者的电脑和手机,发现死者居然都...
    沈念sama阅读 83,932评论 2 376
  • 文/潘晓璐 我一进店门,熙熙楼的掌柜王于贵愁眉苦脸地迎上来,“玉大人,你说我怎么就摊上这事。” “怎么了?”我有些...
    开封第一讲书人阅读 146,770评论 0 330
  • 文/不坏的土叔 我叫张陵,是天一观的道长。 经常有香客问我,道长,这世上最难降的妖魔是什么? 我笑而不...
    开封第一讲书人阅读 53,799评论 1 271
  • 正文 为了忘掉前任,我火速办了婚礼,结果婚礼上,老公的妹妹穿的比我还像新娘。我一直安慰自己,他们只是感情好,可当我...
    茶点故事阅读 62,697评论 5 359
  • 文/花漫 我一把揭开白布。 她就那样静静地躺着,像睡着了一般。 火红的嫁衣衬着肌肤如雪。 梳的纹丝不乱的头发上,一...
    开封第一讲书人阅读 48,069评论 1 276
  • 那天,我揣着相机与录音,去河边找鬼。 笑死,一个胖子当着我的面吹牛,可吹牛的内容都是我干的。 我是一名探鬼主播,决...
    沈念sama阅读 37,535评论 3 390
  • 文/苍兰香墨 我猛地睁开眼,长吁一口气:“原来是场噩梦啊……” “哼!你这毒妇竟也来了?” 一声冷哼从身侧响起,我...
    开封第一讲书人阅读 36,200评论 0 254
  • 序言:老挝万荣一对情侣失踪,失踪者是张志新(化名)和其女友刘颖,没想到半个月后,有当地人在树林里发现了一具尸体,经...
    沈念sama阅读 40,353评论 1 294
  • 正文 独居荒郊野岭守林人离奇死亡,尸身上长有42处带血的脓包…… 初始之章·张勋 以下内容为张勋视角 年9月15日...
    茶点故事阅读 35,290评论 2 317
  • 正文 我和宋清朗相恋三年,在试婚纱的时候发现自己被绿了。 大学时的朋友给我发了我未婚夫和他白月光在一起吃饭的照片。...
    茶点故事阅读 37,331评论 1 329
  • 序言:一个原本活蹦乱跳的男人离奇死亡,死状恐怖,灵堂内的尸体忽然破棺而出,到底是诈尸还是另有隐情,我是刑警宁泽,带...
    沈念sama阅读 33,020评论 3 315
  • 正文 年R本政府宣布,位于F岛的核电站,受9级特大地震影响,放射性物质发生泄漏。R本人自食恶果不足惜,却给世界环境...
    茶点故事阅读 38,610评论 3 303
  • 文/蒙蒙 一、第九天 我趴在偏房一处隐蔽的房顶上张望。 院中可真热闹,春花似锦、人声如沸。这庄子的主人今日做“春日...
    开封第一讲书人阅读 29,694评论 0 19
  • 文/苍兰香墨 我抬头看了看天上的太阳。三九已至,却和暖如春,着一层夹袄步出监牢的瞬间,已是汗流浃背。 一阵脚步声响...
    开封第一讲书人阅读 30,927评论 1 255
  • 我被黑心中介骗来泰国打工, 没想到刚下飞机就差点儿被人妖公主榨干…… 1. 我叫王不留,地道东北人。 一个月前我还...
    沈念sama阅读 42,330评论 2 346
  • 正文 我出身青楼,却偏偏与公主长得像,于是被迫代替她去往敌国和亲。 传闻我的和亲对象是个残疾皇子,可洞房花烛夜当晚...
    茶点故事阅读 41,904评论 2 341

推荐阅读更多精彩内容

  • 引言:马上期末考试了,最近在复习计算机算法分析与程序设计;动态规划,这门课程中最难的几个部分之一,上课老师讲时自己...
    cp_insist阅读 5,174评论 0 3
  • 《算法导论》这门课的老师是黄刘生和张曙,两位都是老人家了,代课很慢很没有激情,不过这一章非常有意思。更多见:iii...
    mmmwhy阅读 5,247评论 5 31
  • 今天来讨论一道算法题,这道算法题我在做的时候真的是没什么思路,甚至看到解法之后依然想了好久才想明白,好久没看过算法...
    柳年思水阅读 14,748评论 0 6
  • 以前我的兴趣爱好太多,这也想做,那也想做,有的是做了没兴趣了,也浪费了时间,现在想来,人这一辈子能把一件事做到最好...
    下雨了收衣服阅读 218评论 2 1
  • 我用一周左右的时间在每天来回半小时的公车上看完了这本书,它讲述的是一个碌碌无为的青年梦想成为一个百万富翁故事。结果...
    巍巍子阅读 1,236评论 0 2