算法竞赛（Python）-状态间的奇妙转移（动态规划）

一、初探动态规划

1 拼图游戏（从搜索到动态规划）

  小余买了一套拼图玩具，玩具包含了n种不同的拼图碎片，每种拼图碎片的宽度都为1，长度各不相同，并且每种拼图碎片的数量足够多（可以认为是无限多）,小余给你出了一个难题，有多少种方式可以拼出1XL的形状。
  输入格式：n表示碎片种类数；L表示要拼出形状的长度。I=[I0,I1,I2…]表示每种碎片的长度。
  输出格式:输出一个整数，表示用n种拼图碎片拼出1XL的方案数。由于数值很大，输出结果保留除以10⁹ +7的余数即可。
  样例输入n=2,L=5;I=[1,2]。样例输出：8。

  假设数据范围比较小，L<=40,从左往右，每一步枚举下一块拼图碎片的形状，就可以用搜索算法遍历每一种可行解。

n=2#碎片种类数
L=5#要拼出形状的长度
I=[1,2]#每种碎片的长度
#用深度优先搜索算法求解length，当前已经拼出的长度
def solve(length):if length>L:return 0elif length==L:return 1 #1种方案else:ans=0#枚举下一块拼图碎片的形状，将其长度加到length中，继续求解for i in range(n):ans+=solve(length+I[i])return ansif __name__=="__main__":print(solve(0))#结果为8

  可惜的是，这种算法通过直接构建可行解来进行计算，所以程序执行的世界与可行解的数量直接相关。如你所见，当把样例中的L改为40时，可行解的数量达到165580141，这样的算法没有办法处理太大的数据。
  既然构造可行解不可行，那么尝试构造“状态”。什么是状态？状态就是通过可行解的中间过程，可以根据实际的问题需要，任意地定义状态。在这个问题中，可以根据已经拼好的形状来定义状态，由于是从左往右拼，所以可以把状态定义为“拼好长度为i的部分”，用1个数组f[i]就可以记录“拼好长度为i部分的方案数”。
  每一个fi如何定义？首先考虑最简单的情况：i=0与i=1,这是两个初始状态，f0=f1=1。
  至于后续的状态，则需要状态之间的转移来计算。当i=2时，要拼好长度为i的部分，考虑最后一块的形状，如果最后一块的形状为1X1，那么就要先拼接好长度为i-1的部分；如果最后一块是1X2,那么就要拼接好长度i-2的部分。所以有：

  脱离这个样式，可以得到一个更一般的公式：   这个方程为状态转移方程，是动态规划的核心所在。

n=2#碎片种类数
L=5#要拼出形状的长度
I=[1,2]#每种碎片的长度
mod=1e9+7
f=[0 for i in range(1005)]
def main():#求解f[0]=1for i in range(1,L+1): #遍历[1,L]#枚举最后一块拼图碎片，利用前面状态的计算结果进行计算for j in range(n):if i-I[j]>=0:f[i]=(f[i]+f[i-I[j]])%modif __name__=="__main__":main()print(f[L])#结果为8

  至此，问题得以解决。动态规划这个词虽然看起来很高端，但实际的代码是非常精简的。
  是不是通过样例得到的公式有点熟悉，这不就是斐波那契数列吗？在前面提到可以用记忆化搜索的方式计算斐波那契数列，那么记忆化搜索也可以解决这个问题。

#用记忆化搜索计算f[i]
def caculate_f(i):#如果已经计算过，直接返回答案if f[i]!=-1:return f[i]else:f[i]=0#枚举最后一块拼图碎片，利用前面状态的计算结果进行计算for j in range(n):#枚举每一种碎片if i-I[j]>=0:f[i]=(f[i]+caculate_f(i-I[j]))%modreturn f[i]
if __name__=="__main__":n = 2  # 碎片种类数L = 5  # 要拼出形状的长度I = [1, 2]  # 每种碎片的长度mod = 1e9 + 7f = [-1 for i in range(1005)] #-1表还未被计算过f[0]=1 #初始化print(caculate_f(L)) #8

  这是动态规划的两种实现方式，不妨把前一种实现方式称为非递归形式，后一种实现方式称为递归形式。这两种实现方式各有千秋，虽然时间复杂度相同，但非递归形式往往快一些。在这个问题中，递归形式显得画蛇添足，但对于一些难以确定状态计算先后顺序的问题，递归形式更容易用代码实现。

2 物流仓库——状态的转移

  通过前面的例子，已初识了动态规划，动态规划的一般步骤如下：

1. 定义状态
2. 确认状态间的转移关系
3. 构造状态转移方程
4. 完成代码

   其中第2步有时会比较困难，接下来这个例子帮你学会厘清状态间的转移关系。
  小余新建了一座物流仓库来临时存放货物。为了节省仓库空间，小余希望把货物堆得高高的，不过货物的形状各不相同，通常把体积大的货物堆放在下层才能保持稳定。已知有n件货物，某些货物可以上下堆放，这些货物要满足以下条件。
  （1）如果货物B可以堆放在货物A上方，那么货物A不能堆放在货物B上方。
  （2）如果货物B可以堆放在货物A上方，同时货物C可以堆放货物B上方，那么货物C可以堆放在货物A上方。
  现在请你帮小余计算这些货物可以堆放多高？
  样例输入格式：n=5;表示货物个数。货物关系[[2,1],[1,4],[3,4],[0,2],[2,4],[0,3]];其中的[a,b]表示货物a可以放在货物b上方。输出4。
  首先定义状态，总共有n个货物，用hi表示货物i放在最上方时最多可以堆放的货物数量。那么要求就是求所有h_i的最大值。
   接下来考虑状态间的转移关系，题目中给出了6个关系，恰好就是状态间的转移关系，如果货物v可以放在货物u上方，那么h_u+1j就是h_v可能的取值，货物v下方堆放的货物显然越多越好，自然可以得到状态转移方程为：

  注意：当计算h_v时，必须保证对每个可以堆放货物v下方的货物u，h_u都已经计算完毕。要保证这样的计算顺序，用记忆化搜索可以实现。

n=5
pre=[[]for i in range(5)]#pre[v]存储所有可以堆放在货物v下方的货物 [[2, 3], [4], [1, 4], [4], []]
pre[2].append(1)
pre[1].append(4)
pre[3].append(4)
pre[0].append(2)
pre[2].append(4)
pre[0].append(3)h=[0 for i in range(10000)] #h[v]表示货物v放在最上方时最多可以堆放的货物数量def max_height(v):if h[v]!=0:return h[v]else:for u in pre[v]:h[v]=max(h[v],max_height(u))h[v]+=1return h[v]if __name__=="__main__":ans=0for i in range(n):ans=max(ans,max_height(i))print(ans) #结果为4

  如果用抽象的眼光看待这个问题，则该问题的本质就是在一个有向图中找到最长的一条链，由于题目中的条件限制，该图不仅仅是一个有向图，还是一个有向无环图，正因为如此，才能顺利用记忆化搜索遍历所有状态。
  在动态规划中，所有的状态关系构成了一个有向无环图，如果状态间的转移关系出现了环，那么就无法使用动态规划来解决问题。
  前面提到，动态规划通常有2种实现方式——递归形式和非递归形式，记忆化搜索就是递归形式，非递归形式也可以解决这个问题，只不过需要计算前将所有状态进行排序，保证在计算每一个货物时，所有可以堆放在其下方的货物都已经计算完毕。下面介绍如何排序。
  开始时，货物4下方不能堆放放在任何货物，没有任何指向货物4，于是直接得到f₄=1。f_i表示i放在最上方时可以堆放的货物数量。

  计算完f₄后,把货物4从图中删除，那么就没有任何边指向货物1和货物3，换句话说，货物1和货物3下方可以堆放的货物都已经计算完了，这时可以计算f₁和f₃。

  或许你已经发现规律了，没有任何边指向的状态，恰是接下来可以计算的状态，每次把这样的状态找出来，计算完毕后将其在图中删除，重复这个过程，就可以找到符合条件的计算顺序。

  这样的排序过程被称为拓扑排序，任何一个有向无环图都可以在符合条件的拓扑排序结果，但拓扑排序的结果不唯一。为了使读者更容易理解拓扑排序的过程，可在代码中把拓扑排序作为一个独立的过程来实现。

n=5
pre=[[]for i in range(5)]#pre[v]存储所有可以堆放在货物v下方的货物 [[2, 3], [4], [1, 4], [4], []]
pre[2].append(1)
pre[1].append(4)
pre[3].append(4)
pre[0].append(2)
pre[2].append(4)
pre[0].append(3)nex=[[]for i in range(5)] #nex[v]存储所有可以堆放在货物v上方的货物[[], [2], [0], [0], [1, 3, 2]]
nex[1].append(2)
nex[4].append(1)
nex[4].append(3)
nex[2].append(0)
nex[4].append(2)
nex[3].append(0)h=[0 for i in range(10000)] #h[v]表示货物v放在最上方时最多可以堆放的货物数量#=============拓扑排序===================
def TopolgicalSort():#in_deg[v]记录了v之前还有多少状态未计算in_deg=[0 for i in range(n)]for u in range(n):for v in nex[u]:in_deg[v]+=1#print(in_deg)[2, 1, 2, 1, 0]#node存储了所有可以计算但未计算的状态node=[]for u in range(n):if in_deg[u]==0:node.append(u)#print(node)#4#permutation存储拓扑排序的结果permutation=[]while len(node)>0:#从node中取出一个可以直接计算的状态u=node[-1] #取出末尾node.pop() #删除最后#将该状态放入排序结果中permutation.append(u)#将该状态从图中删除，更新相关的in_degfor v in nex[u]:in_deg[v]-=1if in_deg[v]==0:node.append(v)return permutation #[4, 3, 1, 2, 0]#===========================动态规划====================
def calculate_max_height(permutation):ans=0for v in permutation:for u in pre[v]:h[v]=max(h[v],h[u])h[v]+=1ans=max(ans,h[v])return ansif __name__=="__main__":permutation=TopolgicalSort()print(calculate_max_height(permutation)) #结果为4

现在，得到了一个重要的结论：动态规划问题中的状态转移关系构成有向无环图，计算时需要安装有向无环图的拓扑排序结果进行。在大多数情况下，动态规划问题的状态关系比较简单，用非递归形式实现较为容易；如果难以厘清状态间的转移关系，可以使用递归形式来实现。

二、状态的巧妙定义

1 不同的状态和转移

  小余想要在股票市场上多赚一些钱，所以不会局限于只进行一次买卖。另一方面，众所周知，过于频繁的买卖操作难以获得稳定的收益，所以，小余想要计算当交易次数（买入次数加卖出次数）不超过q时，能够获得的最大收益是多少。
  具体来说，已知连续n天的股票价格序列为p₀,p₁,…p_n-1 ,开始小余只有1元，假设在第i天小余手里有x元，如何进行买入操作，小余就会获得x/p_i 个单位的股票‘假设在第i天小余手里有y个单位的股票，如果进行卖出操作，小余就会获得yp_i元。买入股票的当天不能卖出，卖出股票的当天不能买入。最后一天结束时，小余手里的资金就是投资的收益，没有卖出的股票不能算作收益。
  输入格式：n=6表示股票交易日的天数，q=4表示最大交易次数；p=[1,2,3,2,1,2]表示股票价格。
  输出格式：输出一个实数，表示小余以1元作为本金，在交易次数不超过q的前提下，可以获得的最大收益。
   为了表示方便，把价格序列p₀,P₁,…P_n-1 的下标修改为从1开始，也就是p₁,p₂，… p_n。
  第i天只有1个状态是不够的，可以将它拆成多个状态，用m(i,j）表示直到第i天为止恰好进行了j次交易能够获得的最大收益。
  n个状态被拆分为nX(q+1）个状态，另外，为了计算方便，添加初始状态m(0,0),m(0,1)…m(0,q)，表示还未进行任何交易时小余手中的资金。公式为：

  这些状态之间的转移关系取决于交易如何进行，考虑利用第k天的状态计算第i天的状态m(i,j),由于定义的状态m(k,j)是指第k天结束时获得的最大收益，所以只能在第k+1天进行交易。

• （1）如果从第k+1天到第i(i>k)天，没有进行任何交易，那么第k天结束时的收益会原封不动地保留到第i天。
• （2）如果在第k+1t天买入，并在第i(i>k)天卖出，那么进行这两次交易后的收益是m(k,j-2)p_i/p_k+1.
• (3)如果在这期间发生了多次买卖，例如，在第k+1天买入，第k+2天卖出，那么就需要利用第k+2天的状态计算m(i,j)，而不是利用第k天的状态计算m(i,j)。

  状态转移方程为：

n=6 #天数
q=4 #最大交易次数
p=[0,1,2,3,2,1,2]#表示股票价格,第0天为0，第1天为1
m=[[0 for i in range(100)]  for i in range(100)] #m[i][j]表示第i天为主恰好进行了j次交易获得的最大收益def max_profit():m[0][0]=1#最开始只有1元for i in range(1,n+1):#遍历每一天[1-6]for j in range(q+1):#交易次数[0-4]for k in range(i):#m[k][j]->m[i][j] 从第k+1天开始到第i天为止不进行任何交易profit=m[k][j]if profit>m[i][j]:m[i][j]=profit#m[k][j-2]->m[i][j] 在第k+1天买入，第i天卖出if (j-2>=0 and k+1!=i):profit=(m[k][j-2]/p[k+1])*p[i]if profit>m[i][j]:m[i][j] = profit#计算最终答案ans=0for j in range(q+1):ans=max(ans,m[n][j])return ansif __name__=="__main__":print(max_profit())#结果6

  不过需要注意，用三层for循环时，时间复杂度达到了O(n² q),但这里需要更低的复杂度，怎么办？答案就是继续拆分状态。
   目前的思路把一次买卖当作一次状态转移，由于需要枚举买卖的时间间隔，才导致了这么高的时间复杂度。现在把每一个状态再次一分为二。

• m(i,j)(Money）:表示直到第i天为止恰好进行了j次交易能够获得的最大收益。
• s(i,j)(Stock):表示直到第i天为止恰好进行了j次交易能够获得的最多股票。

   利用新的状态s(i,j),可以利用新的链式转移关系替代原来的转移关系。

n=6 #天数
q=4 #最大交易次数
p=[0,1,2,3,2,1,2]#表示股票价格,第0天为0，第1天为1
m=[[0 for i in range(100)]  for i in range(100)] #m[i][j]表示第i天为主恰好进行了j次交易获得的最大收益
s=[[0 for i in range(100)]  for i in range(100)] #s[i][j]表示第i天为主恰好进行了j次交易获得的最多股票def max_profit():m[0][0]=1#最开始只有1元for i in range(1,n+1):#遍历天数[1,n]for j in range(q+1):#遍历交易次数#s[i-1][j]->s[i][j] 第i天不买不卖profit=s[i-1][j]if profit>s[i][j]:s[i][j]=profit#s[i-1][j-1]->m[i][j] 第i天卖出if j-1>=0:profit=s[i-1][j-1]*p[i]if profit>m[i][j]:m[i][j]=profit#m[i-1][j-1]->s[i][j] #第i天买入if j-1>=0:profit=m[i-1][j-1]/p[i]if profit>s[i][j]:s[i][j]=profit#m[i-1][j]->m[i][j] 第i天不买不卖profit=m[i-1][j]if profit>m[i][j]:m[i][j]=profit#计算最终答案ans=0for j in range(q+1):ans=max(ans,m[n][j])return ansif __name__=="__main__":print(max_profit()) #结果6

2 流浪猫的家——状态压缩与状态剪枝

   小余和他的朋友收养了一群流浪猫。打算把这群可爱的小猫们安置在后院里。后院可以认为是一片nxn的网格，每个格子可以安置一只小猫。小猫们有很强的领地意识，它们不希望其他的小猫出现在自己附近，所以每个小猫的安置地不能相邻。另外，由于后院中有些杂物，所以有些格子不能安置小猫。那么问题是，后院里最多可以安置多少只流浪猫。

   格式“.”代表该位置可以安置小猫。“#”代码杂草。“c”代表小猫。

   解题思路：把整行放在一起考虑。一行有n个网格，把它们压缩到一起，一共有2ⁿ 种情况，每一种情况按照二进制方式处理可以编码为一个[0,2ⁿ）范围内的整数。下面定义状态。
   f(i,s):第i行采用安置方式s时，前i行可以容纳的最多小猫数量。

   可知共有nx2^n^ 种状态，每一行的2^n^ 状态根据前一行的2^n^ 个状态计算，也就是说，行与行之间有2^n^ x 2^n^ 个转移关系，所有的状态之间就有nx2^n^ x2^n^ 个转移关系，如果n很大，将导致时间复杂度爆炸。   但是，其实结果没有那么糟糕。以n=4为例，看似有16种安置方式，去掉所有存在小猫相邻的安置方式，实际只有8种安置方式。所有合理的安置方式构成集合S，那么状态数总有nxS。

  然后考虑行与行之间的转移关系，由于小猫们不能相邻，所以并非任何两种安置方式都可以放在相邻的两行中。对于每一个安置方式s，可以预先计算出所有可以与它相邻的安置方式集合。有

  与运算，它会一一比较两个二进制数的每一位，如果对应位置上都为1，该位置计算结果也为1，否则为0。利用这种运算符，可以不用将整数分解成二进制形式，而是可以直接快速判断两种安置方式相邻时有没有小猫相邻。

  与运算，还可以用来判断安置方式是否与杂物网格冲突。

import numpy as npn=5 #5行5列
S=[]#不考虑杂物网格时所有合理的安置方式
T=[[] for i in range(np.power(2,n))]#对于每种安置方式可以进行转移的安置方式集合# 判断这一行是否有小猫相邻
def check(state):""":param state: 整数，代表方案:return:"""while(state>0):#示例7的二进制是111,3的二进制是11。7&3==3。7的二进制只取出两位和3进行与运算。因为3只有两位if state & 3==3:return Falsestate=state>>1 #二进制向右移动1位。即删除最右那位二进制。return True
#计算得到所有无小猫相邻的状态
def calculate_S():for i in range(np.int(np.power(2,n))):#遍历方案数if check(i):S.append(i)#计算对于每种安置方式可以进行转移的安置方式集合
def calculate_T():for s in S:for t in S:if s & t ==0:T[s].append(t)if __name__ =="__main__":calculate_S()calculate_T()#输出所有状态转移的次数num=0for s in S:num+=len(T[s])print("状态转移次数：",num*n)#495

  每一种安置方式包含的小猫数量N(s）也可以预先计算出来，N(s）就是整数s用二进制表示的数据。
  状态转移方程已产生，如果安置方式，s与第i行的杂物网格不冲突，那么有:

  否则f(i,s)=0。
  至于如何输出最优安置方式，在计算f(i,s)的过程中，顺便把f(i,s)最大的t记录下来，就可以记录每行每种安置方式前一行应该采取什么样的安置方式。动态规划结束后，用倒序的方式就可以得到所有网格的最优安置方式。

import numpy as npn=5 #5行5列
block=[[0,0,0,0,0],[0,1,0,0,0],[0,1,1,1,0],[1,1,0,1,0],[0,0,0,0,0]] #网格中的杂物情况，1代表有杂物，0代表没有
f=[[0 for i in range(np.power(2,n))] for j in range(n)] #f[i][s]表示第i行采用安置方式s时，前i行可以容纳小猫的最多数量
p=[[0 for i in range(np.power(2,n))] for j in range(n)]#p[i][s]表示f[i][s]取得最大值的前一行安置方式
N=[0 for i in range(np.power(2,n))]#每一种安置方式容纳的小猫数量
S=[]#不考虑杂物网格时所有合理的安置方式
T=[[] for i in range(np.power(2,n))]#对于每种安置方式可以进行转移的安置方式集合# 判断这一行是否有小猫相邻def check(state):""":param state: 整数，代表方案:return:"""while(state>0):#示例7的二进制是111,3的二进制是11。7&3==3。7的二进制只取出两位和3进行与运算。因为3只有两位if state & 3==3:return Falsestate=state>>1 #二进制向右移动1位。即删除最右那位二进制。return True#计算得到所有无小猫相邻的状态
def calculate_S():for i in range(np.int(np.power(2,n))):#遍历方案数if check(i):S.append(i)#计算每种安置方式包含的小猫数量(二进制1的个数] 辅助函数
def count_binary_ones(num):# 此处写你的代码num=bin(num)number=list(num)count=0for i in number:if i=='1':  # 遍历的时候需要注意字符和数字是不同的，需要将字符'1'和数字1比较count+=1return count
#计算每种安置方式包含的小猫数量(二进制1的个数] 主函数
def calculate_N():for s in S:#s为一个整数N[s]=count_binary_ones(s)#计算对于每种安置方式可以进行转移的安置方式集合
def calculate_T():for s in S:for t in S:if s & t ==0:T[s].append(t)#动态规划
def dp():for i in range(n):#将本行的杂物编码为一个整数，注意顺序b=0for j in range(n-1,-1,-1):#[n-1....0]b=b*2+block[i][j]if i==0:#第一行单独考虑for s in S:if s & b==0:f[i][s]=N[s]#f[i][s]表示第i行采用安置方式s时，前i行可以容纳小猫的最多数量.N(s)每一种安置方式容纳的小猫数量else:#从第2行开始，每一行根据上一行的状态进行转移for s in S:if s & b!=0:continue #本行该安置方式与杂物网格冲突for t in T[s]:#T[s]对于每种安置方式可以进行转移的安置方式集合if f[i-1][t]+N[s]>f[i][s]:f[i][s]=f[i-1][t]+N[s]p[i][s]=t#输出
def output():#计算最后一行的最优安置方式state=0for s in S:if f[n-1][s]>f[n-1][state]:state=s#倒序计算每一行的最优安置方式for i in range(n-1,-1,-1):#[n-1....0]#将本行的安置方式填充到block数组中，用2表示小猫temp=statefor j in range(n):if temp & 1:block[i][j]=2temp=temp>>1#二进制向右移动1位。即删除最右那位二进制。state=p[i][state]#输出结果for i in range(n):line = ""for j in range(n):if block[i][j]==0:#没有安置小猫line=line+str(".")+str(" ")elif block[i][j]==1:#杂物line = line + str("#") + str(" ")else:#小猫line = line + str("c") + str(" ")print(line)if __name__ =="__main__":calculate_S()calculate_N()calculate_T()dp()output()print(block)

  在这个问题中，用了两个关键思路：

• 状态压缩：把一行的状态压缩成一个整数，便于用位运算来加速计算。
• 状态剪枝：去掉所有不合理的状态（转移方式），提升算法效率。

三 转移方式的神奇优化

  高效的动态规划算法，不仅需要合理的状态定义，还需要简洁高效的转移方式。

1 运输计划——在转移中剪枝

   一方有难，八方支援。某地突发自然灾害，小余的物流公司要勇于承担起社会责任，尽快将救灾物资运往灾区。
   已知总共有n种物资，每种救援物资分别有m₀, m₁,…m_n-1 箱，每种救援物资每箱的重量分别为w₀，w₁，…w_n-1，每种救援物资每箱的价值分别是v₀,v₁,…v_n-1。小余的物流公司运输能力有限，目前只能将总重量不超过W的物资运往灾区，小余希望最大限度发挥公司的运输能力，在不超重的前提下将总价值尽可能大的物资运往灾区。
   输入格式：n=3表示物资种类。W=7表示最大载重。m=[3,2,2]表示每种物资的箱数。w=[2,3,4]表示每种物资一箱的重量。v=[2,1,2]表示每种物资一箱的价值。
   为了表述方便，把下标修改为从1开始。
   直接解决这个问题可能有点复杂，可以先简化，假设每一种物资都只有一箱，总共有n箱物资，现在要做的就是在这n箱物资中挑选n箱运往灾区。

• 第1步：定义状态。决策的关键因素是物资种类和重量，用f(i,j)表示只考虑前i箱物资时，载重为j时能够运输的最大物资价值。
• 第2步:确定转移关系，.。计算f(i,j)时，需要考虑第i箱物资，如果不运输第i箱物资，那么此时只能在前i-1件物资中挑选，能够运输的最大物资价值就是f(i-1,j)；如果运输第i件物资，那么除了第i件物资以外，还能够运输重量为j-w_i的物资，这部分的价值最大是f(i-1,j-w_i）。所有f(i,j)只需要利用f(i-1,j)和f(i-1,j-w_i)两个状态就可以计算。
• 第3步：构造状态转移方程：
• 

、

n=0 #物资的箱数 初始化
W=0#可以运输的最大重量
w=[0 for i in range(1000)] #每箱物资的重量
v=[0 for i in range(1000)] #每箱物资的价值
f=[[0 for i in range(1000)] for j in range(1000)] #f[i][j]表示只考虑前i箱物资，载重为j时可以运输的最大物资价值def max_value():#依次考虑第i箱物资for i in range(1,n+1):#[1,n]#依次考虑每个可能的重量for j in range(W+1):#[0,W]f[i][j]=f[i-1][j]if j>=w[i]:#当前最大可载重大于物资重量f[i][j]=max(f[i][j],f[i-1][j-w[i]]+v[i])#计算最终的答案ans=0for j in range(W+1):#[0,W]ans=max(ans,f[n][j])return ansif __name__=="__main__":n=int(input("请输入物资箱数n："))W = int(input("请输入可运输最大重量W："))for i in range(1,n+1):#[1,n]w[i]=int(input("请输入第%d箱物资重量w[i]:"%(i)))v[i]=int(input("请输入第%d箱物资价值v[i]:"%(i)))print(max_value())

   时间复杂度为0(nW)，但是这段代码解决的是简化版的问题，要解决原版的问题，一个思路就是把每箱物资单独考虑，只需要在输入数据时进行修改即可。

n=0 #物资的箱数 初始化
W=0#可以运输的最大重量
w=[0 for i in range(1000)] #每箱物资的重量
v=[0 for i in range(1000)] #每箱物资的价值
f=[[0 for i in range(1000)] for j in range(1000)] #f[i][j]表示只考虑前i箱物资，载重为j时可以运输的最大物资价值def max_value():#依次考虑第i箱物资for i in range(1,n+1):#[1,n]#依次考虑每个可能的重量for j in range(W+1):#[0,W]f[i][j]=f[i-1][j]if j>=w[i]:#当前最大可载重大于物资重量f[i][j]=max(f[i][j],f[i-1][j-w[i]]+v[i])#计算最终的答案ans=0for j in range(W+1):#[0,W]ans=max(ans,f[n][j])return ansif __name__=="__main__":package_num=3#物资的种类W=7 #可运输的最大重量m1=[2,3,2]#每种物资的箱数w1=[3,2,4]#每种物资的重量v1=[1,2,2]#每种物资的价值for i in range(package_num):#遍历每种物资for j in range(m1[i]):#遍历当前物资种类的每一箱n += 1  # 箱数加1w[n]=w1[i] #当前物资种类的重量v[n]=v1[i]#当前物资种类的价值print("每箱物资的重量：",w)print("每箱物资的价值：",v)print(max_value())

   只不过这样一来，空间和时间复杂度都爆炸啦。
   先解决空间复杂度，因为前面计算前i箱物资的状态只需用到前i-1箱物资的状态，所以没必要把所有的状态都存下来，只需前i箱物资的状态，动态地更新，用新的状态覆盖原来的状态，则用一个一维数组即可。

f=[0 for i in range(1000)] #f[j]表示载重为j时可以运输的最大物资价值def max_value():#依次考虑第i箱物资for i in range(1,n+1):#依次考虑每个可能的载重量for j in range(W,-1,-1):#[W,0]if j>=w[i]:#载重量大于当前箱物资重量f[j]=max(f[j],f[j-w[i]]+v[i])#计算最终的答案ans=0for j in range(W+1):#ans=max(ans,f[j])return ans

   接下来考虑时间的问题，如果一种物资有x箱，那么就需要考虑运输0箱、1箱、…x箱，这x+1种情况需要计算x+1次，这是导致时间复杂度爆炸的原因。
   可以将这些物资打包运算，以15箱为例，利用二进制的原理，把它们打包成4个包裹，分布包含1箱，2箱，4箱，8箱，这4包物资就可以组合成0,1…15中的每一个数字。对4包物资进行计算远比15箱物资进行计算快很多。更一般情况，把x箱物资打包包含2⁰ ,2¹ ,2² …箱包裹，最后剩下不足2ⁿ的部分单独打包。
   全部代码

n=0 #物资的箱数 初始化
W=0#可以运输的最大重量
w=[0 for i in range(1000)] #每箱物资的重量
v=[0 for i in range(1000)] #每箱物资的价值
f=[0 for i in range(1000)] #f[j]表示载重为j时可以运输的最大物资价值def max_value():#依次考虑第i箱物资for i in range(1,n+1):#依次考虑每个可能的载重量for j in range(W,-1,-1):#[W,0]if j>=w[i]:#载重量大于当前箱物资重量f[j]=max(f[j],f[j-w[i]]+v[i])#计算最终的答案ans=0for j in range(W+1):#ans=max(ans,f[j])return ansif __name__=="__main__":package_num=3#物资的种类W=7 #可运输的最大重量m1=[2,3,2]#每种物资的箱数w1=[3,2,4]#每种物资的重量v1=[1,2,2]#每种物资的价值for i in range(package_num):#遍历每种物资#以二进制形式拆分物资k=1while k<=m1[i]:m1[i]-=kn += 1  # 箱数加1w[n] = k*w1[i]  # 当前物资种类的重量v[n] = k*v1[i]  # 当前物资种类的价值k=k*2#如果还有剩余，单独打包剩余物资if m1[i]>0:w[n] = m1[i]*w1[i]  # 剩余这种类物资的所有重量v[n] = m1[i]*v1[i]  # 剩余这种类物资的所有价值print(max_value())#结果6

2 会议安排——在决策中剪枝

   小余的物流公司上市了，为了确定公司未来的发展规划，小余作为董事长，决定召开一次会议。包括小余在内的n位公司成员参与这次会议。
   关于公司未来的发展规划，每个人都有自己的想法，如果n个人参与会议，那么任何两个人之间都免不了有一场争论，必须等所有人争论完之后会议才能结束。举例说明，如果有3个人参会，那么0与1,1与2,0与2之间发生三场争论，耗时分别为t_0,1 ,t_1,2 ,t_0,2 ,整个会议需要t_0,1 +t_1,2 +t_0,2才能结束。
   为了节省时间，小余决定让会议分组进行，把n个人分成k组，每组进行一场会议，依次进行所有的k场会议，总共所需要的时间是这k场会议所需的时间之和，需要注意的是，只能将连续相邻的几位公司成员分到同一组。请你帮小余确定分组，尽快结束会议。
   输入格式
n=8 #参与会议的人数
m=4#分组数
t=[[0,0,0,1,1,1,1,1],
[0,0,0,1,1,1,1,1],
[0,0,0,1,1,1,1,1],
[1,1,1,0,1,1,1,1],
[1,1,1,1,0,1,1,1],
[1,1,1,1,1,0,2,2],
[1,1,1,1,1,2,0,2],
[1,1,1,1,1,2,2,0]]#t[i,j]表示成员i与成员j被分到同一组时争论的时间

   在这个问题中，有两个影响决策的关键因素——分组数和人数，转态就按照分组数和人数定义，用f(i,j)表示把前j个人（1,2…j）分到前i组中需要的最短会议时间。
   至于状态转移，则是遍历第i组的分组情况，利用前i-1组的状态计算前i组的状态，例如，当前把前k个会议成员分到前i-1组中时，把会议成员k+1,k+2,…j分配到第i组，就可以把前j个会议成员分到前i组了，换句话说，从状态f(i,k)转移到了状态f(i,j)。也就有了状态转移方程。

   其中的T(k+1,j）表示将会议成员k+1,k+2,…j分到第i组进行会议所需的时间，其实就是矩阵{u_i,j}中子矩阵和的一半。如图所示：【第k+1行到第j行，第k+1列到第j列的数据之和为2倍T(k+1,j）因为每个数字被计算了2次】

   为了快速计算T(k+1,j），定义s(i,j)为矩阵{u_i,j}中第i行第j列的元素及其左上方所有元素之和，利用{是（s(i,j)}可以快速巧妙地计算T(k+1,j），如图所示：

n=8 #参与会议的人数
m=4#分组数
t=[[0,0,0,1,1,1,1,1],[0,0,0,1,1,1,1,1],[0,0,0,1,1,1,1,1],[1,1,1,0,1,1,1,1],[1,1,1,1,0,1,1,1],[1,1,1,1,1,0,2,2],[1,1,1,1,1,2,0,2],[1,1,1,1,1,2,2,0]]#t[i,j]表示成员i与成员j被分到同一组时争论的时间 序号从0开始#为了表述方便，修改t[i,j]坐标从1开始。
u=[[0 for i in range(n+2)] for j in range(n+2)]#u[i,j]表示成员i与成员j被分到同一组时争论的时间 序号从1开始
for i in range(1,n+1):for j in range(1,n+1):u[i][j]=t[i-1][j-1]s=[[0 for i in range(n+2)] for j in range(n+2)]#u矩阵的前缀和
f=[[0 for i in range(n+2)] for j in range(n+2)]#表示把前j个人分到前i组中需要的最短会议时间#计算成员l，l+1,...r分到同一组进行会议需要的时间
def T(l,r):"""就是公式中的T(k+1,j):param l: int:param r: int:return:"""return (s[r][r]-s[l-1][r]-s[r][l-1]+s[l-1][l-1])/2
# 动态规划
def dp():#把f[i][j]初始化为足够大的数据inf=1e9for i in range(m+1):#[0,m]遍历各个分组for j in range(n+1):#遍历各个人f[i][j]=inf#初始状态，0个人分到0个组开会时间为0分钟f[0][0]=0#在状态转移中计算f[i][j]for i in range(1,m+1):#[1,m]遍历各个分组for j in range(i,n+1):#遍历各个分组，1个分组至少1人for k in range(i-1,j):#[i-1,j-1]if f[i-1][k]+T(k+1,j)<f[i][j]:f[i][j]=f[i-1][k]+T(k+1,j)if __name__=="__main__":#计算前缀和sfor i in range(1,n+1):#遍历每个人for j in range(1,n+1):s[i][j]=u[i][j]+s[i-1][j]+s[i][j-1]-s[i-1][j-1]#动态规划dp()#输出print(f[m][n]) #结果为3

  观察上面的代码，时间复杂度为O(mn²）,最核心也是最耗时的部分是dp()函数内部的三层for循环，最内存的for循环是在穷举每一个可能的k。找到最优的分组方案。这一步的本质是极其“暴力”的穷举算法，下面尝试缩小穷举的范围。
  在动态规划中，把每一个状态f(i,j)最优的k记录下来，记作d(i,j),如果最优的k有多个，只记录最小的那一个，用数学公式表达出来就是：

  可以上图中发现规律，每一行都是单调递增的，每一列也都是单调递增的。即：

  这个结论其实不难理解，一组中包含的人数越多，对应的子矩阵范围越大，分组会议花费的时间越多。为了尽快结束所有会议，分组时要尽可能分得均匀一点。所以，当分组个数或人数增加时，最优的划分点有向右移动的趋势。
  在穷举每一个k时，可以直接把范围[i-1,j-1]缩小到[d(i-1,j),d(i,j+1)]。最内层循环执行的次数是≤mn+n² 。
  需要注意的是，由于计算d(i,j)需要用到d(i-1,j)与d(i,j+1)，所以需要调整计算顺序，第2层循环中的j要从大往小遍历。

n=8 #参与会议的人数
m=4#分组数
t=[[0,0,0,1,1,1,1,1],[0,0,0,1,1,1,1,1],[0,0,0,1,1,1,1,1],[1,1,1,0,1,1,1,1],[1,1,1,1,0,1,1,1],[1,1,1,1,1,0,2,2],[1,1,1,1,1,2,0,2],[1,1,1,1,1,2,2,0]]#t[i,j]表示成员i与成员j被分到同一组时争论的时间 序号从0开始#为了表述方便，修改t[i,j]坐标从1开始。
u=[[0 for i in range(n+2)] for j in range(n+2)]#u[i,j]表示成员i与成员j被分到同一组时争论的时间 序号从1开始
for i in range(1,n+1):for j in range(1,n+1):u[i][j]=t[i-1][j-1]s=[[0 for i in range(n+2)] for j in range(n+2)]#u矩阵的前缀和
f=[[0 for i in range(n+2)] for j in range(n+2)]#表示把前j个人分到前i组中需要的最短会议时间d=[[0 for i in range(n+2)] for j in range(n+2)]#f[i][j]取得最优解的k#计算成员l，l+1,...r分到同一组进行会议需要的时间
def T(l,r):"""就是公式中的T(k+1,j):param l: int:param r: int:return:"""return (s[r][r]-s[l-1][r]-s[r][l-1]+s[l-1][l-1])/2
# 动态规划
def dp():#把f[i][j]初始化为足够大的数据inf=1e9for i in range(m+1):#[0,m]遍历各个分组for j in range(n+1):#遍历各个人f[i][j]=inf#初始化d[i][j]的边界for i in range(1,n+1):#[1,n]d[0][i]=0d[i][n+1]=n#初始状态，0个人分到0个组开会需要0分钟f[0][0]=0#在状态转移中计算f[i][j]for i in range(1,m+1):#遍历每个分组for j in range(n,i-1,-1):#[n,i]for k in range(int(max(i-1,d[i-1][j])),int(min(j-1,d[i][j+1]))+1):#[max(i-1,d[i-1][j])),min(j-1,d[i][j+1]]if f[i-1][k]+T(k+1,j)<f[i][j]:f[i][j]=f[i-1][k]+T(k+1,j)d[i][j]=kif __name__=="__main__":#计算前缀和sfor i in range(1,n+1):#遍历每个人for j in range(1,n+1):s[i][j]=u[i][j]+s[i-1][j]+s[i][j-1]-s[i-1][j-1]#动态规划dp()#输出print(f[m][n]) #结果为3

三、经典的动态规划算法

  严格来讲，动态规划是一种思想。而非一种算法，把用动态规划思想构造出的算法统称为动态规划算法。

1 路径规划——用动态规划创造算法

  小余的物流公司旗下有n个物流网点，编号为0,1,…n-1，这些物流网点之间有若干条道路连接，有些道路是双向通行的，有些道路是单向通行的，小余急需将一批货物从物流网点s运输到哦物流网点t，需规划出最短运输路径。

  输入格式:第1行是三个整数n,m₁,m₂，分别表示物流网点的个数、双向通行道路的条数、单向通行道路的条数;第二行是两个整数s,t，表示运输的起点和终点;接下来m₁行，每行三个整数u,v,w，表示物流网点u与v之间有一条长度为w的双向通行道路;接下来m₂行，每行三个整数u,v,w，表示从物流网点u到v有一条长度为w的单向通行道路。
  输出格式:
  输出一个整数，表示从物流网点s到t的最短路径长度。

样例输入：
4 1 4
1 2
0 1 1
1 3 4
0 3 2
3 2 3
2 0 4
样例输出：6

  首先，有3个显然成立但非常有用的结论：

• 1 所有道路都可以被认为是单向通行的。双向通行的道路可以被认为是两条单向通行的道路。
• 2 最短路径上不可能绕圈。如果从物流网点s到t的某一条路径上经过某个位置u至少两次，那么这条路径一定不是s到t的最短路径。
• 3 只有最短路径才能构成新的最短路径。如果从物流网点s到t的最短路径上经过某个位置u，那么可以把这段路径拆分成s->u,u->t两部分，这两部分一定是s到u，u到t的最短路径。

  动态规划

• 第1步：定义状态f(i,j,k)表示从物流网点i到j，中间只经过前k个点的最短路径长度。最初k=0，如果没有从点i到点j的道路，那么f(i,j,0)=+∞（或者一个足够大的数据），否则f(i,j,0)就是点i到点j的道路中最短的长度。
  

• 第2步：确定转移关系，从只经过前k-1点，到只经过前k个点，就是向现有的最短路径中添加第k个点。对应任何两个点i,j，有两种路径可能成为新的最短路径，一种是不走点k，最短路径的长度是f(i,j,k-1);另一种是走点k，那么就要先到点k，再从点k到点j，只有最短路径才能构成新的最短路径，所以两短路径的长度分别是f(i,k,k-1)与f(k,j,k-1)。
  

• 第3步：构造状态转移方程：
  

   时间复杂度O(n³）,空间复杂度O(n³）、经过前k个点的状态只与经过前k-1个点的状态有关，所以不必用三维数组来存储所有状态，用一个二维数组即可。

n=4 #物流网点个数
m1=1#双向通行道路的个数
m2=4#单向通行道路的个数
s=1;t=2#s起点，t终点
f=[[0 for i in range(10)] for j in range(10)] #f[i][j]表示从点i到j的最短路径长度
#初始化
for i in range(n):for j in range(n):#如果点i到j没有路径，可以认为点i到j有一条很长的路径f[i][j]=1e9f[i][i]=0
#输入所有的双向通行路径f[0][1]=min(f[0][1],1)
f[1][0] = min(f[1][0], 1)
#输入所有的单向通行路径
f[1][3]=min(f[1][3],4)
f[0][3]=min(f[0][3],2)
f[3][2]=min(f[3][2],3)
f[2][0]=min(f[2][0],4)#======================
def Floyd():#依次考虑第k个点for k in range(n):for i in range(n):for j in range(n):#如果把点k插入(i->j)的最短路径中能找到一条更短的路径#那就更新f[i][j]f[i][j]=min(f[i][j],f[i][k]+f[k][j])#输出
Floyd()
print(f[s][t]) #结果6

  这就是著名的Floyd算法，一个用动态规划思想巧妙构造出的、极其简洁的算法。Floyd算法可以直接计算出任何两个点之间的最短路径长度，这既是优点也是缺点，要想计算从点s到点t的最短路径，必须繁琐地把任何两个点之间的最短路径算出来，下面尝试换一种计算方法。

• 第1步：定义状态。用d(j,k)表示从起点s出发，至多经过k条路径到达点j的最短路径长度，最开始k=0，只有d(s,0)=0,其余的d(j,0)=+∞。从起点s出发，到任何一个点的最短路径上不会出现超过n-1条路径，否则就会出现绕圈的情况。如图所示：最终点s到点t的路径就是d(t,n-0)
  

• 第2步，确定转移关系。从至多经过k-1条路径，到至多经过k条路径，就是向现有的最短路径末尾添加一条路径。对于状态d(j,k),如果选择用原来的路径，那么长度为d(j,k-1);如果选择在某个路径末尾添加一条长度为w的路径u->j，先到点u再到点j，那么长度为d(u,k-1)+w。

• 第3步：构造状态转移方程，即：式子中w(u->j)是路径u->j的长度。
  

n=4 #物流网点个数
m1=1#双向通行道路的个数
m2=4#单向通行道路的个数
s=1;t=2#s起点，t终点
m=2*m1+m2 #道路的个数
u=[0 for i in range(20)] #每条道路的起点
v=[0 for i in range(20)] #每条道路的终点
w=[0 for i in range(20)] #每条道路的长度
d=[0 for i in range(20)]#d[i]表示从起点s出发到达点i的路径长度# 输入所有的双向通行路径
u[0]=0;v[0]=1;w[0]=1#第1条路径[0,1]=1
u[1] = 1;v[1] = 0;w[1] = 1  # 第2条路径[1,0]=1
#输入所有的单向通行路径
u[2] = 1;v[2] = 3;w[2] = 4  # 第3条路径[1,3]=4
u[3] = 0;v[3] = 3;w[3] = 2  # 第4条路径[0,3]=2
u[4] = 3;v[4] = 2;w[4] = 3  # 第5条路径[3,2]=3
u[5] = 2;v[5] = 0;w[5] = 4  # 第6条路径[2,0]=4def Bellman_Ford():#初始化inf=1e9for i in range(n):#遍历物流网点d[i]=infd[s]=0 #起点#重复n-1次for k in range(n-1):#枚举每一条道路，将其衔接到某条最短路径末尾for i in range(m):d[v[i]]=min(d[v[i]],d[u[i]]+w[i])#求解
Bellman_Ford()
print(d[t])#结果6

   这是另一个著名的Bellman-Ford算法，如果只需求两点之间的最短路径。在路径比较稀疏时，比Floyd算法快得多。

2 矩阵乘积——用动态规划优化算法

   众所周知，矩阵乘法是线性代数中的基本运算，在机器学习等领域有非常广泛的应用。矩阵A与矩阵B相乘时，需要保证矩阵A的列数等于矩阵B的行数。假设A是n行k列的矩阵，B是k行m列的矩阵，那么AB是n行m列的矩阵，其中第i行第j列元素恰好是矩阵A的第i行与矩阵B的第j列对应位置相乘再相加的结果，所以，朴素的矩阵乘法计算的时间复杂度为O(nkm)，其中n,k,m分别是矩阵A的行数，矩阵A的列数和矩阵B的列数。
   为了适当简化问题，这里认为矩阵A与矩阵B相乘所需要的计算次数恰好是nkm。假设为了满足某些计算需求，需要计算一串矩阵的乘积A₀A₁A₂…A_n-1,其中第i个矩阵的行数和列数分别是m_i和m_i+1,从左往右依次相乘，需要的计算次数是:m₀m₁m₂+m₀m₂m₃+m₀m₃m₄+…+m₀m_n-1m_n。
   矩阵乘法虽然不满足交换律，但满足结合律。如A₀*A₁*A₂*A₄ ，可以写成A₀ *（A₁*A₂）*A₄

  输入格式：第1行是一个整数n，表示矩阵的个数；第2行是n+1个整数m₀,m₁,m₂,…m_n,其中第i个矩阵的行数和列数分别是m_i和m_i+1。输出格式：输出一个整数，表示计算这n个矩阵的乘积需要的最少计算次数。
样例输入:n=4;m=[3,4,5,2,3].输出82。
  调整矩阵乘法的计算顺序，无非就是在其中添加括号，一个括号就是一段区间，很容易想到用区间定义状态。用f(l,r)表示计算A_lA_l+1…A_r需要的最少次数。
  以上矩阵乘积的计算可以认为是区间合并的过程。如图所示，最初有n个长度为1的区间，每次挑选相邻的两个区间合并，经过n-1次区间合并后，矩阵的乘积计算完毕。

  最初的计算顺序就是最优的区间合并顺序，对于某一个区间[l,r]来说，通过合并得到区间[l.r]之间，必然是两个区间[l,i]与[i+1，r]。要算出区间[l,r]对应的矩阵乘积，需要的最少的计算次数恰好是f(l.r);要算出[i+1，r]对应的矩阵乘积，需要的最少的计算次数是f(i+1,r)。这两段区间的计算结果分别是m_l行m_i+1列，m_i+1行m_r+1列的矩阵。计算这两个矩阵的乘积需要m_lm_i+1m_r+1次计算。所以，转态转态方程是：

n=4#矩阵个数
m=[3,4,5,2,3]#矩阵的行数和列数
f=[[0 for i in range(10)] for j in range(10)]#f[l][r]表示计算A[l]A[l+1]...A[r]需要的最少计算次数
def dp():#注意for循环中的计算顺序for r in range(n):for l in range(r-1,-1,-1):#[r-1,0]f[l][r]=1e9for i in range(l,r):f[l][r]=min(f[l][r],f[l][i]+f[i+1][r]+m[l]*m[i+1]*m[r+1])#求解
dp()
#输出
print(f[0][n-1])#结果82