-
2019-07-11 16:16:24
动态规划(Dynamic programming)
是一种在数学、计算机科学和经济学中使用的,通过把原问题分解为相对简单的子问题的方式求解复杂问题的方法。 动态规划算法是通过拆分问题,定义问题状态和状态之间的关系,使得问题能够以递推(或者说分治)的方式去解决。
什么是动态规划
动态规划(Dynamic Programming)对于子问题重叠的情况特别有效,因为它将子问题的解保存在表格中,当需要某个子问题的解时,直接取值即可,从而避免重复计算!
动态规划是一种灵活的方法,不存在一种万能的动态规划算法可以解决各类最优化问题(每种算法都有它的缺陷)。所以除了要对基本概念和方法正确理解外,必须具体问题具体分析处理,用灵活的方法建立数学模型,用创造性的技巧去求解。
基本策略
基本思想与分治法类似,也是将待求解的问题分解为若干个子问题(阶段),按顺序求解子阶段,前一子问题的解,为后一子问题的求解提供了有用的信息。在求解任一子问题时,列出各种可能的局部解,通过决策保留那些有可能达到最优的局部解,丢弃其他局部解。依次解决各子问题,最后一个子问题就是初始问题的解。
动态规划中的子问题往往不是相互独立的(即子问题重叠)。在求解的过程中,许多子问题的解被反复地使用。为了避免重复计算,动态规划算法采用了填表来保存子问题解的方法。
适用问题
那么什么样的问题适合用动态规划的方法来解决呢?
适合用动态规划来解决的问题,都具有下面三个特点:最优化原理、无后效性、有重叠子问题。
(1)最优化原理:如果问题的最优解所包含的子问题的解也是最优的,就称该问题具有最优子结构,即满足最优化原理。
(2)无后效性:即某阶段状态一旦确定,就不受这个状态以后决策的影响。也就是说,某状态以后的过程不会影响以前的状态,只与当前状态有关。
(3)有重叠子问题:即子问题之间是不独立的,一个子问题在下一阶段决策中可能被多次使用到。(该性质并不是动态规划适用的必要条件,但是如果没有这条性质,动态规划算法同其他算法相比就不具备优势。
这类问题的求解步骤通常如下:
初始状态→│决策1│→│决策2│→…→│决策n│→结束状态
(1)划分:按照问题的特征,把问题分为若干阶段。注意:划分后的阶段一定是有序的或者可排序的
(2)确定状态和状态变量:将问题发展到各个阶段时所处的各种不同的客观情况表现出来。状态的选择要满足无后续性
(3)确定决策并写出状态转移方程:状态转移就是根据上一阶段的决策和状态来导出本阶段的状态。根据相邻两个阶段状态之间的联系来确定决策方法和状态转移方程
(4)边界条件:状态转移方程是一个递推式,因此需要找到递推终止的条件
算法实现
动态规划三要素:
(1)问题的阶段
(2)每个阶段的状态
(3)相邻两个阶段之间的递推关系
整个求解过程可以用一张最优决策表来描述,最优决策表是一张二维表(行:决策阶段,列:问题的状态)表格需要填写的数据一般对应此问题的在某个阶段某个状态下的最优值(如最短路径,最长公共子序列,最大价值等),填表的过程就是根据递推关系,最后根据整个表格的数据通过简单的取舍或者运算求得问题的最优解。
例如:f(n,m)=max{f(n-1,m), f(n-1,m-w[n])+P(n,m)}
1 事例一背包问题
问题描述:假设我们有n种类型的物品,分别编号为1, 2...n。其中编号为i的物品价值为vi,它的重量为wi。为了简化问题,假定价值和重量都是整数值。现在,假设我们有一个背包,它能够承载的重量是Cap。现在,我们希望往包里装这些物品,使得包里装的物品价值最大化,那么我们该如何来选择装的东西呢?注意:每种物品只有一件,可以选择放或者不放。初始化数据为:n=5,w={2,2,6,5,4},v={6,3,5,4,6},Cap=10
解法如下:
A)描述最优解的结构
设子问题:f[i][v]表示允许前i件物品放入容量为v的背包时可以获得的最大价值。注:这里的i从0到5,v从0到10
为了能够得到已经计算过的,更小规模的子问题,我们可以根据当前限重来只考虑第i件物品放或者不放,那么就可以转化为涉及前i-1件物品的问题,
即:
情况1、如果第i件物品不能放(即这个物品的重量直接大于了当前限重v),则问题转化为“前i-1件物品放入容量为v的背包中”,即f[i-1][v];
情况2、如果放第i件物品是可以放也可以不放,则问题转化为:
1)、如果选择不放第i件物品,则问题转化为“前i-1件物品放入容量为v的背包中”,即变大时f[i-1][v];
2)、如果选择放第i件物品,则问题转化为“前i-1件物品放入剩下的容量为v-w[i]的背包中”,此时能获得的最大价值就是f [i-1][v-w[i]]再加上通过放入第i件物品获得的价值v[i]。
则情况2下,f[i][v]的值就是1),2)中最大的那个值。
最优子结构描述如下:当子问题f[i][v]是最优时,其子问题f[i-1][v]和f[i-1][v-w[i]](中的较大者)显然同样也必须是最优的值,不然在情况1或者情况2下总会矛盾。
B) 递归定义最优解的值
根据上面的分析,显然子问题
f[i][v]=f[i-1][v],这时是情况1
f[i][v]=max{f[i-1][v], f[i-1][v-w[i]]+v[i] },这时是情况2。
C)按自底而上的方式计算最优解的值
import numpy as np #行李数n,不超过的重量W,重量列表w和价值列表p def fun(n,W,w,p): a=np.array([[0]*(W+1)]*(n+1)) #依次计算前i个行李的最大价值,n+1在n的基础上进行 for i in range(1,n+1): for j in range(1,W+1): if w[i-1]>j: a[i,j]=a[i-1,j] else: a[i,j]=max(a[i-1,j],p[i-1]+a[i-1,j-w[i-1]])#2种情况取最大值 #print(a) print('max value is'+str(a[n,W])) findDetail(p,n,a[n,W]) #找到价值列表中的一个子集,使得其和等于前面求出的最大价值,即为选择方案 def findDetail(p,n,v): a=np.array([[True]*(v+1)]*(n+1)) for i in range(0,n+1): a[i][0]=True for i in range(1,v+1): a[0][i]=False for i in range(1,n+1): for j in range(1,v+1): if p[i-1]>j: a[i,j]=a[i-1,j] else: a[i,j]=a[i-1,j] or a[i-1,j-p[i-1]] if a[n,v]: i=n result=[] while i>=0: if a[i,v] and not a[i-1,v]: result.append(p[i-1]) v-=p[i-1] if v==0: break i-=1 print(result) else: print('error') weights=[1,2,5,6,7,9] price=[1,6,18,22,28,36] fun(len(weights),13,weights,price)
2 事例二
有n级台阶,一个人每次上一级或者两级,问有多少种走完n级台阶的方法。
分析:动态规划的实现的关键在于能不能准确合理的用动态规划表来抽象出 实际问题。在这个问题上,我们让f(n)表示走上n级台阶的方法数。
那么当n为1时,f(n) = 1,n为2时,f(n) =2,就是说当台阶只有一级的时候,方法数是一种,台阶有两级的时候,方法数为2。那么当我们要走上n级台阶,必然是从n-1级台阶迈一步或者是从n-2级台阶迈两步,所以到达n级台阶的方法数必然是到达n-1级台阶的方法数加上到达n-2级台阶的方法数之和。即f(n) = f(n-1)+f(n-2),我们用dp[n]来表示动态规划表,dp[i],i>0,i<=n,表示到达i级台阶的方法数。
class Solution: """ @param n: an integer @return: an ineger f(n) """ def up(self, n): # write your code here # if n == 0: # return 0 L = [] L.append(1) L.append(2) for i in range(2, n): L.append(L[i - 1] + L[i - 2]) return L[n - 1]
3 事例三
给定一个矩阵m,从左上角开始每次只能向右走或者向下走,最后达到右下角的位置,路径中所有数字累加起来就是路径和,返回所有路径的最小路径和,如果给定的m如下,那么路径1,3,1,0,6,1,0就是最小路径和,返回12.
1 3 5 9
8 1 3 4
5 0 6 1
8 8 4 0
分析:对于这个题目,假设m是m行n列的矩阵,那么我们用dp[m][n]来抽象这个问题,dp[i][j]表示的是从原点到i,j位置的最短路径和。我们首先计算第一行和第一列,直接累加即可,那么对于其他位置,要么是从它左边的位置达到,要么是从上边的位置达到,我们取左边和上边的较小值,然后加上当前的路径值,就是达到当前点的最短路径。然后从左到右,从上到下依次计算即可。
m[i][j] = min(m[i-1][j],m[i][j-1]) + p[i][j]
代码不再详细写,与背包问题类似的矩阵求解。
4 事例四
给定两个字符串str1和str2,返回两个字符串的最长公共子序列,例如:str1="1A2C3D4B56",str2="B1D23CA45B6A","123456"和"12C4B6"都是最长公共子序列,返回哪一个都行。
分析:本题是非常经典的动态规划问题,假设str1的长度为M,str2的长度为N,则生成M*N的二维数组dp,dp[i][j]的含义是str1[0..i]与str2[0..j]的最长公共子序列的长度。
dp值的求法如下:
dp[i][j]的值必然和dp[i-1][j],dp[i][j-1],dp[i-1][j-1]相关,结合下面的代码来看,我们实际上是从第1行和第1列开始计算的,而把第0行和第0列都初始化为0,这是为了后面的取最大值在代码实现上的方便,dp[i][j]取三者之间的最大值。
# 此例中有多个相同长度的公共子串,但只能获取第一个子串 def find_lcsubstr(s1, s2): # 下面4行不要直接写在循环中,减少计算 s1_len = len(s1) + 1 #为方便后续计算,多了1行1列 s2_len = len(s2) + 1 s3_len = len(s1) s4_len = len(s2) m = [[0 for j in range(s2_len)] for i in range(s1_len)] #生成0矩阵 maxNum = 0 #初始最长匹配长度 p = 0 #匹配的起始位置 for i in range(s3_len): for j in range(s4_len): if s1[i] == s2[j]: #相同则累加 m[i + 1][j + 1] = m[i][j] + 1 #给相同字符赋值,值为左上角值加1 if m[i + 1][j + 1] > maxNum: maxNum = m[i + 1][j + 1] #获取最大匹配长度 p = i + 1 #记录最大匹配长度的终止位置 print(m) return s1[p - maxNum : p], maxNum #返回最长子串及其长度 print(find_lcsubstr(str_a, str_b))
5 事例五
把只包含质因子2、3和5的数称作丑数(Ugly Number)。例如6、8都是丑数,但14不是,因为它包含质因子7。 习惯上我们把1当做是第一个丑数。求按从小到大的顺序的第N个丑数。
解题思路:这个用动态规划感觉就很好啦
终止条件:第N个丑数
状态:s[i] 第i个丑数
状态转移方程: s[i] = min(s[t1]*2 , s[t2]*3, s[t3]*5)。这个t1,t2,t3<i-1
这3个t怎么确定呢?从i-1往后倒腾,满足s[t1-1]*2小于等于s[i-1]但s[t1]*2大于s[i-1]
确定了之后敲代码把~~(代码还有很多可以优化的,但觉得这个可读性比较强,贴上自己的代码)
class Solution: def GetUglyNumber_Solution(self, index): if index<1: return 0 if index==1: return 1 s = [] s.append(1) t1 = 0 t2 = 0 t3 = 0 for i in range(1,index): for j in range(i): if s[j]*2 > s[i-1]: t1 = j break for j in range(i): if s[j]*3 > s[i-1]: t2 = j break for j in range(i): if s[j]*5 > s[i-1]: t3 = j break s.append(min(s[t1]*2,s[t2]*3,s[t3]*5)) return s[-1]
6 事例六
找零钱问题,已经零钱面额为1、3、4,求找零n所用零钱数最少的方案
解答过程:对于找零n的最少零钱数f(n),它和f(n-1),f(n-3),f(n-4)有关,即它等于这3者中最小的值加1.
# 找零钱字典,key为面额,value为最少硬币数 change_dict = {} def rec_change(M, coins): change_dict[0] = 0 s = 0 for money in range(1, M+1): num_of_coins = float('inf') #意思是要求50的最少找零数,在46,47,49的最少找零数中找到最小的即可 for coin in coins: if money >= coin: # 记录每次所用的硬币数量 if change_dict[money-coin]+1 < num_of_coins: num_of_coins = change_dict[money-coin]+1 s = coin #记录每次找零的面额 change_dict[money] = num_of_coins return change_dict[M],s # 求出具体的找零方式 # 用path变量记录每次找零的面额 def method(M, coins): print('Total denomination is %d.'%M) nums, path = rec_change(M, coins)#path为最少硬币数方案中的一个面额值 print('The smallest number of coins is %d.'%nums) print('%s'%path, end='') while M-path > 0: M -= path nums, path = rec_change(M, coins) print(' -> %s'%path, end='') print() coins = (1, 3, 4) method(50, coins)
7 事例七
钢条切割,已经各长度的钢条和对应的收益,问长度为n的钢条怎么切割收益最大。
要求长度为n的钢条切割最大收益,则在n-1最大收益+长度1的收益,n-2最大收益+长度2最大收益……中取最大者。那么依次求长度1~n的钢条最大收益即可。
# 钢条长度与对应的收益 length = (1, 2, 3, 4,5, 6, 7, 8, 9, 10) profit = (1, 5, 8, 9,10, 17, 17, 20, 24, 30) # 参数:profit: 收益列表, n: 钢条总长度 def bottom_up_cut_rod(profit, n): r = [0] # 收益列表 s = [0]*(n+1) # 切割方案列表 for j in range(1, n+1): q = float('-inf') #每次循环求出长度为j的钢条切割最大收益r[j],s[j]则保存切割方案中最长的那一段长度 for i in range(1, j+1): if max(q, profit[length.index(i)]+r[j-i]) == profit[length.index(i)]+r[j-i]:#元组index从1开始 s[j] = i#如果切割方案为1和2,那么2会覆盖1,即保存最长的一段 q = max(q, profit[length.index(i)]+r[j-i]) r.append(q) #r[n]保存长度为n钢条最大切割收益 return r[n], s[n] # 切割方案 def rod_cut_method(profit, n): how = [] while n != 0: t,s = bottom_up_cut_rod(profit, n) how.append(s) n -= s return how #输出长度1~10钢条最大收益和最佳切割方案 for i in range(1, 11): t1 = time.time() money,s = bottom_up_cut_rod(profit, i) how = rod_cut_method(profit, i) t2 = time.time() print('profit of %d is %d. Cost time is %ss.'%(i, money, t2-t1)) print('Cut rod method:%s\n'%how)
8 事例八
水杯摔碎问题,有n个水杯和k层楼,求最少测试几次可以确定水杯刚好在哪一层楼摔碎。
解答过程:假设从x层楼开始扔为f(n,x),如果水杯碎了水杯数量-1需要探测的楼层为x-1层,则为f(n-1,x-1),如果没碎水杯还是n个需要探测k-x层,则为f(n,k-x)
import numpy as np #n个水杯k层楼,最少需要几次测试确定水杯在几层楼刚好摔破 def solvepuzzle(n, k): numdrops = np.array([[0]*(k+1)]*(n+1)) for i in range(k+1): numdrops[1, i] = i#只有一个水杯,最坏情况是跟楼层数一样 for i in range(2, n+1):#2到n个水杯 for j in range(1, k+1):#楼层1到k minimum = float('inf') #每次循环得出一种(i,j)下的最少次数 for x in range(1, j+1): minimum = min(minimum, (1+max(numdrops[i, j-x], numdrops[i-1, x-1]))) numdrops[i, j] = minimum print(numdrops) return numdrops[n,k] t = solvepuzzle(3, 10) print(t)
9 事例九
给定n个水杯和d次尝试机会,求最多能探测多少楼层。
解答过程:f(d,n)=f(d-1,n)+f(d-1,n-1),令g(d,n)=f(d,n+1)-f(d,n)=g(d-1,n)+g(d-1,n-1),这跟二项式C(n,k)=C(n-1,k)+C(n-1,k-1)相似,故g(d,n)=C(d,n)-->f(d,n)=求和(C(d,i)) i从1到n-1,i>=d时C(d,i)=0
#n个水杯d次尝试机会,最多探测多少层楼? #f(d,n)=求和i=1~n-1{C(d,i)} 对所有d>=1 and i<d def assist(d,i):#C(d,i) sum=1 sum2=1 for k in range(d-i+1,d+1): sum*=k for j in range(1,i+1): sum2*=j sum=sum/sum2 return sum def f(d,n): if d<1 or n<1: return 0 elif d==1: return 1 sum=0 for i in range(1,n+1): if i<d: sum+=assist(d,i) return int(sum) print(f(3,3)) #这里可以逆向求解n个水杯k层楼至少需要多少次测试 def reverseFun(n,k): d=1 while f(d,n)<k: d+=1 print(d) return d reverseFun(3,10)
更多相关内容 -
运筹系列8:Set partitioning问题的列生成方法
2018-08-15 22:32:47列生成算法是不断添加求解变量的算法,可参考论文。列生成算法常常用于如下的场景:使用set-covering构建的模型,变量非常多,约束相对较少。 具体来说,场景如下:有nnn个0-1变量z1...znz1...znz_1...z_n,每个...1. 应用场景
列生成算法是不断添加求解变量的算法,可参考论文。上一节我们介绍了约束条件存在分块结构的情形下可以使用单纯型法的判断规则进行列生成,下面我们介绍另一种场景,set-partitionning问题的列生成算法。
具体来说,场景如下:有 n n n个0-1变量 z 1 . . . z n z_1...z_n z1...zn,每个 z i z_i zi带着很多中间变量 x i , j x_{i,j} xi,j用来进行约束,这是个变量很多,约束也很多的模型。我们首先对问题模型进行转化,将所有 z i z_i zi的组合枚举出来,在枚举的过程中,就把约束条件都过了一遍,只有满足所有约束条件的组合才会保留下来。最后的枚举结果分别对应 y 1 . . . y m y_1...y_m y1...ym。
举个VRP问题,每一列是一个可能路径(列生成指的就是生成闭环路径的规则),主问题的约束条件则是除了depot外的每个点访问过一次。
再举个航空排班的例子。假设机组有x个人需要执飞y个航班,每一列是一个可能的部分航班计划,而主问题的约束条件则是每个航班执飞过一次。不难看出, m m m是 n n n的指数形式,使用set-partitioning模型之后,变量的数量大的可怕。好消息是,原问题的约束条件比较紧,只有少数组合成立,可以通过一定的规则产生,这样m的数量比较小,问题就比较好解决。
2. 算法步骤
问题的列生成算法对应的模型是:
min Σ j y j \min \Sigma_{j} y_j minΣjyj
s.t. A y ≥ b T Ay ≥ b^T Ay≥bT
A A A是一个矩阵, b = [ b 1 , . . . , b m ] b=[b_1,...,b_m] b=[b1,...,bm]和 y = [ y 1 , . . . , y n ] y=[y_1,...,y_n] y=[y1,...,yn]是向量;并且有重要的一点:set-covering的变量 y j y_j yj对应的列 a j = [ a 1 j , . . . , a m j ] T a_j=[a_{1j},...,a_{mj}]^T aj=[a1j,...,amj]T生成的规则,可以通过线性约束条件 B ∗ a j ≤ D B*a_j ≤ D B∗aj≤D得到。这使得我们可以用线性规划的方法来求解子问题~
注意由于列生成算法求解的是松弛线性规划问题,因此对整数规划模型需要结合branch and bound算法进行。2.1 限制主问题
首先用启发式算法找出一部分 A A A,比如说选出了 k k k列。然后我们的线性规划问题就变成了:
min y 1 + y 2 + . . . + y k y_1+y_2+...+y_k y1+y2+...+yk
s.t. a i , 1 ∗ y 1 + a i , 2 ∗ y 2 + . . . + a i , k ∗ y k ≥ b i a_{i,1}*y_1+a_{i,2}*y_2+...+a_{i,k}*y_k≥ b_i ai,1∗y1+ai,2∗y2+...+ai,k∗yk≥bi, i = 1... n i = 1...n i=1...n
相比原来的模型,相当于把 y k + 1 y_{k+1} yk+1~ y m y_m ym强制限制为非基变量了,称为限制主问题(Restricted Master Problem,RMP)。上面的限制主问题求解完成后,我们想使用单纯型法进行基变量的转换,看看 y k + 1 y_{k+1} yk+1~ y m y_m ym中,是否有可以转入基变量的列。检验数 σ = c N − π ∗ a j \sigma = c_N - \pi * a_j σ=cN−π∗aj,并且 π = c B B − 1 π = c_BB^{-1} π=cBB−1可以由求解器给出。我们要找出非基变量中最小的负数 σ \sigma σ,将其转入基变量。正如前面所述,我们这里使用线性规划来找这个新的列。2.2 子问题
注意 c = [ 1 , . . . , 1 ] c = [1,...,1] c=[1,...,1],因此 min c N − π a j c_N - \pi a_j cN−πaj等价于
min 1 − π ∗ a 1 - \pi *a 1−π∗a
s.t. B ∗ a ≤ D B*a ≤ D B∗a≤D (列生成的规则)
称为子问题(sub problem)。如果目标函数最优值<0,就将新生成的列yk+1转入基变量,生成新的限制主问题进行求解。如此往复,直至子问题的目标函数≥0。3. 例子
3.1 问题描述
来看网上流传很多的cutting stock problem的例子,问题如下:
我们有一些纸筒,每个纸筒16m。顾客需要裁剪的长度为:25个3m,20个6m,15个7m。要求在满足顾客需求的情况下,裁剪的纸筒数最小。3.2 数学模型和初始可行解
我们可以用启发式算法找到一个upper bound的初始解:
5个方案1:5个3m
10个方案2:2个6m
8个方案3:2个7m
总计23筒。也就是说,我们用23筒是肯定可以满足要求的,这算问题的一个上界。下面我们探索一下其他的切筒方式,看能不能给出下界。
用set-covering对问题进行建模:设 P P P是所有可行的裁剪方案的集合(其中前3个可以设置为我们前面的3个裁剪方案),里面方案的总数为 n n n(我们并不需要确切的知道这个值是多少,只需要知道它很大)。令 a i j a_{ij} aij表示第 j j j种方案里类别 i i i的个数, y j y_j yj表示第 j j j种方案的选择个数,原问题可以变为:
min y 1 + . . . + y n \min y_1+...+y_n miny1+...+yn
a 1 , 1 y 1 + a 1 , 2 y 2 + . . . + a 1 , n y n ≥ 25 a_{1,1}y_1+a_{1,2}y_2+...+a_{1,n}y_n \geq 25 a1,1y1+a1,2y2+...+a1,nyn≥25
a 2 , 1 y 1 + a 2 , 2 y 2 + . . . + a 2 , n y n ≥ 20 a_{2,1}y_1+a_{2,2}y_2+...+a_{2,n}y_n \geq 20 a2,1y1+a2,2y2+...+a2,nyn≥20
a 3 , 1 y 1 + a 3 , 2 y 2 + . . . + a 3 , n y n ≥ 15 a_{3,1}y_1+a_{3,2}y_2+...+a_{3,n}y_n \geq 15 a3,1y1+a3,2y2+...+a3,nyn≥15
其中
a 1 , 1 = 5 a_{1,1}=5 a1,1=5, a 2 , 1 = 0 a_{2,1}=0 a2,1=0, a 3 , 1 = 0 a_{3,1}=0 a3,1=0
a 1 , 2 = 0 a_{1,2}=0 a1,2=0, a 2 , 2 = 2 a_{2,2}=2 a2,2=2, a 3 , 2 = 0 a_{3,2}=0 a3,2=0
a 1 , 3 = 0 a_{1,3}=0 a1,3=0, a 2 , 3 = 0 a_{2,3}=0 a2,3=0, a 3 , 3 = 2 a_{3,3}=2 a3,3=2
对应最初的3种方案。
a i , j a_{i,j} ai,j满足条件: 3 ∗ a 1 , j + 6 ∗ a 2 , j + 7 ∗ a 3 , j ≤ 16 3*a_{1,j}+6*a_{2,j}+7*a_{3,j}\leq 16 3∗a1,j+6∗a2,j+7∗a3,j≤16(前面所说的列生成的规则)
注意我们在使用列生成算法的时候,将整数变量松弛为了连续变量。3.3 第一轮
问题对应的初始解为 y 1 = 5 y_1=5 y1=5, y 2 = 10 y_2=10 y2=10, y 3 = 8 y_3=8 y3=8, y 4 = . . . y n = 0 y_4=...y_n=0 y4=...yn=0
限制主问题如下:
min y 1 + y 2 + y 3 \min y_1+y_2+y_3 miny1+y2+y3
5 y 1 + 0 y 2 + 0 y 3 ≥ 25 5 y_1+0 y_2+0y_3 \geq 25 5y1+0y2+0y3≥25
0 y 1 + 2 y 2 + 0 y 3 ≥ 20 0y_1+2y_2+0y_3 \geq 20 0y1+2y2+0y3≥20
0 y 1 + 0 y 2 + 2 y 3 ≥ 15 0y_1+0y_2+2y_3 \geq 15 0y1+0y2+2y3≥15下面是pymprog的sensitivity输出:
对应 π = ( 0.2 , 0.5 , 0.5 ) \pi = (0.2, 0.5, 0.5) π=(0.2,0.5,0.5),我们把新转入的变量记作 a 4 = ( a 1 , 4 , a 2 , 4 , a 3 , 4 ) a_4 = (a_{1,4},a_{2,4},a_{3,4}) a4=(a1,4,a2,4,a3,4),则子问题为:
min 1 − 0.2 ∗ a 1 , 4 − 0.5 ∗ a 2 , 4 − 0.5 ∗ a 3 , 4 1 - 0.2*a_{1,4}-0.5*a_{2,4}-0.5*a_{3,4} 1−0.2∗a1,4−0.5∗a2,4−0.5∗a3,4
s.t. 3 ∗ a 1 , 4 + 6 ∗ a 2 , 4 + 7 ∗ a 3 , 4 ≤ 16 3*a_{1,4}+6*a_{2,4}+7*a_{3,4}\leq 16 3∗a1,4+6∗a2,4+7∗a3,4≤16 (列生成的规则)
a i , j ∈ Z a_{i,j}\in Z ai,j∈Z
求解结果为 a 4 = ( 1 , 2 , 0 ) a_4=(1,2,0) a4=(1,2,0),检验数 σ = 1 − π ∗ a 4 = − 0.2 ≤ 0 \sigma = 1 - \pi * a_4 =-0.2\le 0 σ=1−π∗a4=−0.2≤0。此时问题目标为22.5。
将 y 4 y_4 y4进基,把这个结果添加到主问题当中去,开始第二轮迭代。3.5 第二轮
添加了 y 4 y_4 y4后,限制主问题为:
min y 1 + y 2 + y 3 + y 4 \min y_1+y_2+y_3+y_4 miny1+y2+y3+y4
5 y 1 + 0 y 2 + 0 y 3 + 1 y 4 ≥ 25 5 y_1+0 y_2+0y_3+1y_4 \geq 25 5y1+0y2+0y3+1y4≥25
0 y 1 + 2 y 2 + 0 y 3 + 2 y 4 ≥ 20 0y_1+2y_2+0y_3+2y_4 \geq 20 0y1+2y2+0y3+2y4≥20
0 y 1 + 0 y 2 + 2 y 3 + 0 y 4 ≥ 15 0y_1+0y_2+2y_3+0y_4 \geq 15 0y1+0y2+2y3+0y4≥15直观上来看,是添加了一个新列,所以称为列生成算法。上面的主问题求解结果为:
对应 π = ( 0.2 , 0.4 , 0.5 ) \pi = (0.2, 0.4, 0.5) π=(0.2,0.4,0.5),我们把新转入的变量记作 a 5 = ( a 1 , 5 , a 2 , 5 , a 3 , 5 ) a_5 = (a_{1,5},a_{2,5},a_{3,5}) a5=(a1,5,a2,5,a3,5),则子问题为:
min 1 − 0.2 ∗ a 1 , 5 − 0.4 ∗ a 2 , 5 − 0.5 ∗ a 3 , 5 1 - 0.2*a_{1,5}-0.4*a_{2,5}-0.5*a_{3,5} 1−0.2∗a1,5−0.4∗a2,5−0.5∗a3,5
s.t. 3 ∗ a 1 , 5 + 6 ∗ a 2 , 5 + 7 ∗ a 3 , 5 ≤ 16 3*a_{1,5}+6*a_{2,5}+7*a_{3,5}\leq 16 3∗a1,5+6∗a2,5+7∗a3,5≤16 (列生成的规则)
a i , j ∈ Z a_{i,j}\in Z ai,j∈Z
求解结果为 a 5 = ( 3 , 0 , 1 ) a_5=(3,0,1) a5=(3,0,1),检验数 σ = 1 − π ∗ a 5 = − 0.1 < 0 \sigma = 1 - \pi * a_5=-0.1< 0 σ=1−π∗a5=−0.1<0。此时问题目标为20.5。
将 y 5 y_5 y5进基,把这个结果添加到主问题当中去,开始第三轮迭代。3.6 第三轮
添加了 y 5 y_5 y5后,限制主问题为:
min y 1 + y 2 + y 3 + y 4 + y 5 \min y_1+y_2+y_3+y_4+y_5 miny1+y2+y3+y4+y5
5 y 1 + 0 y 2 + 0 y 3 + 1 y 4 + 3 y 5 ≥ 25 5 y_1+0 y_2+0y_3+1y_4 +3y_5 \geq 25 5y1+0y2+0y3+1y4+3y5≥25
0 y 1 + 2 y 2 + 0 y 3 + 2 y 4 + 0 y 5 ≥ 20 0y_1+2y_2+0y_3+2y_4 +0y_5\geq 20 0y1+2y2+0y3+2y4+0y5≥20
0 y 1 + 0 y 2 + 2 y 3 + 0 y 4 + 1 y 5 ≥ 15 0y_1+0y_2+2y_3+0y_4 +1y_5\geq 15 0y1+0y2+2y3+0y4+1y5≥15上面的主问题求解结果为:
对应 π = ( 0.166667 , 0.416667 , 0.5 ) \pi = (0.166667, 0.416667, 0.5) π=(0.166667,0.416667,0.5),我们把新转入的变量记作 a 6 = ( a 1 , 6 , a 2 , 6 , a 3 , 6 ) a_6 = (a_{1,6},a_{2,6},a_{3,6}) a6=(a1,6,a2,6,a3,6),则子问题为:
min 1 − 0.166667 ∗ a 1 , 6 − 0.416667 ∗ a 2 , 6 − 0.5 ∗ a 3 , 6 1 - 0.166667*a_{1,6}-0.416667*a_{2,6}-0.5*a_{3,6} 1−0.166667∗a1,6−0.416667∗a2,6−0.5∗a3,6
s.t. 3 ∗ a 1 , 6 + 6 ∗ a 2 , 6 + 7 ∗ a 3 , 6 ≤ 16 3*a_{1,6}+6*a_{2,6}+7*a_{3,6}\leq 16 3∗a1,6+6∗a2,6+7∗a3,6≤16 (列生成的规则)
a i , j ∈ Z a_{i,j}\in Z ai,j∈Z
求解结果为 a 6 = ( 1 , 1 , 1 ) a_6=(1,1,1) a6=(1,1,1),检验数 σ = 1 − π ∗ a 6 = − 0.083334 < 0 \sigma = 1 - \pi * a_6=-0.083334< 0 σ=1−π∗a6=−0.083334<0。此时问题目标为20。
将 y 6 y_6 y6进基,把这个结果添加到主问题当中去,开始第三轮迭代。3.7 第四轮
添加了 y 6 y_6 y6后,限制主问题为:
min y 1 + y 2 + y 3 + y 4 + y 5 + y 6 \min y_1+y_2+y_3+y_4+y_5+y_6 miny1+y2+y3+y4+y5+y6
5 y 1 + 0 y 2 + 0 y 3 + 1 y 4 + 3 y 5 + y 6 ≥ 25 5 y_1+0 y_2+0y_3+1y_4 +3y_5+y_6 \geq 25 5y1+0y2+0y3+1y4+3y5+y6≥25
0 y 1 + 2 y 2 + 0 y 3 + 2 y 4 + 0 y 5 + y 6 ≥ 20 0y_1+2y_2+0y_3+2y_4 +0y_5+y_6\geq 20 0y1+2y2+0y3+2y4+0y5+y6≥20
0 y 1 + 0 y 2 + 2 y 3 + 0 y 4 + 1 y 5 + y 6 ≥ 15 0y_1+0y_2+2y_3+0y_4 +1y_5+y_6\geq 15 0y1+0y2+2y3+0y4+1y5+y6≥15上面的主问题求解结果为:
对应 π = ( 0.2 , 0.4 , 0.4 ) \pi = (0.2, 0.4, 0.4) π=(0.2,0.4,0.4),我们把新转入的变量记作 a 7 = ( a 1 , 7 , a 2 , 7 , a 3 , 7 ) a_7 = (a_{1,7},a_{2,7},a_{3,7}) a7=(a1,7,a2,7,a3,7),则子问题为:
min 1 − 0.2 ∗ a 1 , 7 − 0.4 ∗ a 2 , 7 − 0.4 ∗ a 3 , 7 1 - 0.2*a_{1,7}-0.4*a_{2,7}-0.4*a_{3,7} 1−0.2∗a1,7−0.4∗a2,7−0.4∗a3,7
s.t. 3 ∗ a 1 , 7 + 6 ∗ a 2 , 7 + 7 ∗ a 3 , 7 ≤ 16 3*a_{1,7}+6*a_{2,7}+7*a_{3,7}\leq 16 3∗a1,7+6∗a2,7+7∗a3,7≤16 (列生成的规则)
a i , j ∈ Z a_{i,j}\in Z ai,j∈Z
求解结果为 a 7 = ( 3 , 1 , 0 ) a_7=(3,1,0) a7=(3,1,0),检验数 σ = 1 − π ∗ a 6 = 0 \sigma = 1 - \pi * a_6= 0 σ=1−π∗a6=0,结束迭代,此时问题目标为19。3.8 结果
使用下列代码求出最优解:
from pymprog import * p = model('example') y = p.var('y',6,int) # variables p.minimize(sum(y), 'master') r1 = 5*y[0]+ y[3]+3*y[4]+y[5]>= 25 r2 = 2*y[1]+2*y[3]+y[5] >= 20 r3 = 2*y[2] +y[4]+y[5]>= 15 p.solve() for i in range(6): print(y[i].primal)
求解结果为: y 1 = 1 , y 4 = 3 , y 5 = 1 , y 6 = 14 y_1=1,y_4=3,y_5=1,y_6=14 y1=1,y4=3,y5=1,y6=14,总计需要19个纸筒。
也可以将最后一轮的 y y y取整得到一个可行解: y 1 = 2 , y 4 = 3 , y 6 = 15 y_1=2,y_4=3,y_6=15 y1=2,y4=3,y6=15,需要20个纸筒。
这里给出第四轮的限制主问题和子问题的代码供参考:from pymprog import * p = model('master problem') y = p.var('y',6) # variables p.minimize(sum(y), 'master') r1 = 5*y[0]+ y[3]+3*y[4]+y[5]>= 25 r2 = 2*y[1]+2*y[3]+y[5] >= 20 r3 = 2*y[2] +y[4]+y[5]>= 15 p.solve() p.sensitivity()
from pymprog import * s = model('sub problem') a = s.var('a',3) # variables for i in range(3): a[i].kind = int s.minimize(1-(0.2*a[0] + 0.4*a[1] +0.4*a[2]), 'sub') 3*a[0] +6*a[1] +7*a[2]<= 16 # mountain bike limit s.solve() print(s.vobj()) for i in range(3): print(a[i].primal)
然后下面给出直接求最优解的代码。pymprog用了一天都没有求出来,因此改用cplex进行求解:
import cplex from cplex.exceptions import CplexError my_obj = [0.0,0.0,0.0,1.0]*23 my_ub = [5.0,2.0,2.0,1.0]*23 my_lb = [0.0,0.0,0.0,0.0]*23 my_ctype = "IIII"*23 my_colnames = [] my_rhs = [25,20,15]+[0,16]*23 my_sense = "GGG"+"LL"*23 for i in range(23): my_colnames = my_colnames + ["x0,"+str(i),"x1,"+str(i),"x2,"+str(i),"y"+str(i)] def populatebyrow(prob): prob.objective.set_sense(prob.objective.sense.minimize) prob.variables.add(obj=my_obj, lb=my_lb, ub=my_ub, types=my_ctype,names=my_colnames) row0 = [1.0,0.0,0.0,0.0]*23 row1 = [0.0,1.0,0.0,0.0]*23 row2 = [0.0,0.0,1.0,0.0]*23 rows = [[my_colnames, row0], [my_colnames, row1], [my_colnames, row2]] for i in range(23): row3 = [0.0,0.0,0.0,0.0]*23 row3[i*4:i*4+4]=[1.0,1.0,1.0,-100.0] row4 = [0.0,0.0,0.0,0.0]*23 row4[i*4:i*4+4]=[3.0,6.0,7.0,0.0] rows = rows+[[my_colnames, row3], [my_colnames, row4]] prob.linear_constraints.add(lin_expr=rows, senses=my_sense,rhs=my_rhs) try: my_prob = cplex.Cplex() handle = populatebyrow(my_prob) my_prob.solve() except CplexError as exc: print(exc) print("Solution status = ", my_prob.solution.status[my_prob.solution.get_status()]) print("Solution value = ", my_prob.solution.get_objective_value()) x = my_prob.solution.get_values() for i in range(23): if x[i*4+3]>0: print(x[i*4:i*4+3])
结果如下:
CPXPARAM_Read_DataCheck 1 Tried aggregator 1 time. MIP Presolve modified 23 coefficients. Reduced MIP has 49 rows, 92 columns, and 230 nonzeros. Reduced MIP has 23 binaries, 69 generals, 0 SOSs, and 0 indicators. Presolve time = 0.00 sec. (0.12 ticks) Found incumbent of value 22.000000 after 0.00 sec. (0.48 ticks) Probing time = 0.00 sec. (0.01 ticks) Tried aggregator 1 time. Reduced MIP has 49 rows, 92 columns, and 230 nonzeros. Reduced MIP has 23 binaries, 69 generals, 0 SOSs, and 0 indicators. Presolve time = 0.00 sec. (0.21 ticks) Probing time = 0.00 sec. (0.01 ticks) MIP emphasis: balance optimality and feasibility. MIP search method: dynamic search. Parallel mode: deterministic, using up to 4 threads. Root relaxation solution time = 0.00 sec. (0.13 ticks) Nodes Cuts/ Node Left Objective IInf Best Integer Best Bound ItCnt Gap * 0+ 0 22.0000 0.0000 100.00% 0 0 6.6667 38 22.0000 6.6667 30 69.70% 0 0 10.2791 39 22.0000 Cuts: 62 143 53.28% 0 0 18.0972 40 22.0000 Cuts: 80 213 17.74% * 0+ 0 21.0000 18.0972 13.82% 0 0 18.4688 35 21.0000 Cuts: 38 305 12.05% 0 0 18.7500 19 21.0000 Cuts: 41 395 10.71% 0 0 18.7500 12 21.0000 Cuts: 18 444 10.71% 0 0 18.8667 10 21.0000 Cuts: 12 496 10.16% * 0+ 0 19.0000 18.8667 0.70% 0 0 cutoff 19.0000 18.8667 496 0.70% Elapsed time = 0.18 sec. (7.56 ticks, tree = 0.01 MB, solutions = 3) Implied bound cuts applied: 22 Mixed integer rounding cuts applied: 35 Zero-half cuts applied: 5 Gomory fractional cuts applied: 7 Root node processing (before b&c): Real time = 0.19 sec. (7.57 ticks) Parallel b&c, 4 threads: Real time = 0.00 sec. (0.00 ticks) Sync time (average) = 0.00 sec. Wait time (average) = 0.00 sec. ------------ Total (root+branch&cut) = 0.19 sec. (7.57 ticks) Solution status = MIP_optimal Solution value = 19.0 [1.0, 1.0, 1.0] [1.0, 1.0, 1.0] [1.0, 2.0, 0.0] [1.0, 1.0, 1.0] [1.0, 1.0, 1.0] [1.0, 1.0, 1.0] [1.0, 1.0, 1.0] [1.0, 1.0, 1.0] [1.0, 1.0, 1.0] [1.0, 1.0, 1.0] [1.0, 2.0, 0.0] [1.0, 1.0, 1.0] [1.0, 2.0, 0.0] [1.0, 1.0, 1.0] [5.0, 0.0, 0.0] [1.0, 1.0, 1.0] [1.0, 1.0, 1.0] [3.0, 0.0, 1.0] [1.0, 1.0, 1.0]
最终是19个。
-
算法学习(动态规划)- 数塔问题
2020-11-28 17:58:21想想或许是鸡蛋问题对我现在而言难了一些,所以只好找了一些其它的动态规划问题,从简单入手,先去理解“动态规划”的思想精髓,再反过来去思考“鸡蛋问题”。 其中一个经典问题是01背包问题,这是我之前一直想搞懂...前言
之前碰到了扔鸡蛋问题(给你2个鸡蛋,在100层楼上扔,要求想出一个好的策略,去测出哪一层楼开始鸡蛋就会碎掉),一直摸不着头脑。后来才知道可以使用“动态规划”这种思想(或者叫算法范式)去解决这个问题。
但是看了一些鸡蛋问题和动态规划的文章,依然只是流于形式,并不能理解其中的精髓。想想或许是鸡蛋问题对我现在而言难了一些,所以只好找了一些其它的动态规划问题,从简单入手,先去理解“动态规划”的思想精髓,再反过来去思考“鸡蛋问题”。
其中一个经典问题是01背包问题,这是我之前一直想搞懂的一个问题。看了一篇文章,就大致上理解了这个问题的解决办法和思路。所以我就不班门弄斧了,直接推荐这篇文章:动态规划之01背包问题(最易理解的讲解)。
比起一开始就搬概念和公式,这篇文章用了一个“填表格”的方式去让读者理解背包问题的解决过程,这对初学者更友好。推荐大家也试着去画出那个表格,然后自己填一填,不需要编程,只需要一张纸和笔和你的脑子,就能按照文中的思路把表格填完,然后就解决了背包问题。算法编程实际只是对你脑子的这种思维过程的一个模拟。
我想,之所以这种所谓“思想”、“范式”之类的东西难以理解,也许是因为它经过了双重的“抽象”吧。比如,以下是动态规划思想的描述:
动态规划算法与分治法类似,其基本思想也是将待求解问题分解成若干个子问题,先求解子问题,然后从这些子问题的解得到原问题的解。与分治法不同的是,适合于用动态规划求解的问题,经分解得到子问题往往不是互相独立的。若用分治法来解这类问题,则分解得到的子问题数目太多,有些子问题被重复计算了很多次。如果我们能够保存已解决的子问题的答案,而在需要时再找出已求得的答案,这样就可以避免大量的重复计算,节省时间。我们可以用一个表来记录所有已解的子问题的答案。不管该子问题以后是否被用到,只要它被计算过,就将其结果填入表中。这就是动态规划法的基本思路。具体的动态规划算法多种多样,但它们具有相同的填表格式…
有时候,真的想吐槽这些东西就是“正确的废话”。你不能说这些描述有什么错,但是真遇到问题的时候,光看这些描述你又完全没办法想出解决问题的方法。
解决某个具体问题,我们会把问题进行抽象,例如把搜索问题抽象成树问题。而这些所谓的“思想”,为了归纳这一类问题的共性,又成了对这类问题的抽象。即,成了“抽象之抽象”,双重抽象。
普适的代价即抽象
不知道别人如何,反正我对于“抽象”这类反人类大脑结构的东西实在无法理解。只能重新“化抽象为形象”,才能去“识”得它。就算是概念、思想这类看不见摸不着的东西,还是得转换成我脑海一个具象化的东西才行。
理解了01背包问题之后,我开始找到第2个动态规划的问题,即数塔问题。企图通过解决这些问题,一窥“动态规划”的冰山一角。
下面进入正题。
问题
如图,有一个数塔,从顶部出发,在每一个节点可以选择向左或者向右走,一直走到底层。要求找出一条路径,使得所走路径上的节点数字之和最大。
(图片来自百度)分析
1.穷举
1)思路来源
这是在没有去看正确的动态规划方法之前,我自己想出的思路。
虽然不是“动态规划”,但是这个走歪的方向,也许更能让我们理解究竟怎样才算是“动态规划”。所以这个思路也是有价值的,可以说一说。
我一开始的理解,“动态规划”就是“使用已经计算好的信息,避免重复的计算”。这个想法的来源来自于这篇文章:
文章用了一张图来解释什么是动态规划:
我用当时的个人理解来复述一下:要计算A -> Q -> Z这条线路的最优选择,实际上并不需要我们把 A -> H -> Q -> Z、 A -> I -> Q -> Z 这两条路径都计算出来。 我们在计算A -> Q 的时候, 实际上已经把A -> H ->Q 、A -> I -> Q 这部分计算过了,因此已经得出了A -> Q的最优线路。计算A -> Q -> Z的时候,只需要利用这个结果即可。A -> Q是A -> Q -> Z的一个部分(注意:这个说法不完全正确,所以才会导致我得出了一个不完美的思路),因此能够“使用已经计算好的信息,避免重复的计算”。当时我的理解不能说十分准确,认为动态规划就是:一个子问题如果属于一个大问题的一部分,就可以利用子问题的结果信息去计算大问题,避免重复计算。这个说法不准确,我后面会纠正这种想法。
但此时我们就顺着这种错误想法,看看这种想法会产生怎样的解题思路。
2)思维过程
毫无疑问,我们可以把问题图中的数塔看成一个大数塔,它其实可以看作包含了很多小数塔。比如只有两层的 9、12、15:
或者三层的9、12、15、10、6、8:
总之一个大数塔是可以拆分成小数塔的。那么,为了找到一个大数塔哪条路最优,是不是可以对小数塔找到最优?然后每一层每一层我们都最优,最后就能得到大数塔的最优了呢?咋听起来很诱人,似乎找到了很厉害很巧妙的办法。但可惜,在这里行不通。这实际上就是“贪心算法”(em…又一种算法思想范式…)的本质:只要每一步我们都找最优,找到最后我们就能找到最后的最优。
还是引用推荐文章的一张图:
去计算A -> Z的最优线路,我们先一步一步来。先算A到“H or I”的最优线路,然后往下算“H or I”到“P or Q”的最优线路,最后算“P or Q”到Z的最优线路。算到最后,我们就得到了A到Z的最优线路。这种情况下,我们才能用“贪心法”。看上去似乎跟“动态规划”的图不太一样,我把这张图改改,把“H or I”、“P or Q”这两个节点展开:
可以看到,任意一层的任意节点,都与下一层的所有节点有连接。第一层的A,与第二层的H、I都是连接的。第二层的H、I,与三层的P、Q都是连接的。第三层的P、Q,与第四层的Z都是连接的。大概只有问题的推理形式具有这样的图特性,我们才能使用“贪心”法去解决吧。
我来再来看看我们的数塔,如果我把数塔横着倒放,变成这样:
看起来是不是跟“贪心法”的结构图不一样?是不是更像“动态规划”的那张结构图:
这也就说明了,为什么数塔问题可以用“动态规划”去解决。3)解决办法(当然,最后想歪了,得出了穷举法。。。)
既然没办法通过一步步计算最优,来得出最后的最优,那我就妥协一些。每一次计算某一级子树塔的时候,我就用低一层级的子数塔的全部路线,然后再跟这一层的全部节点作一个连接,得出新的路线。
这样一步步不断计算出路线、以及这些路线的节点数值之和,最后把所有路线全部计算出来,作一个排序,就能知道哪条路线节点数值之和最大,不是就能得出问题答案了么?
想出这个办法的时候,我并不能确定这到底是不是“动态规划”,毕竟这时候对“动态规划”的理解还太浅了。但是又感觉这似乎不是“穷举”,因为我并没有把从根节点到末尾节点的路线全部列举一遍(至少从程序上我不会),我是“一步一步”计算的路线,并且我也把大数塔拆成来小数塔解决,相当于“把大问题拆分成相似的子问题”。从这些特征来看,这似乎也不算是“太笨的方法”,怎么会是“穷举”呢?( ̄▽ ̄)"
再根据这个思路写写程序,发现我居然可以用递归!?这么一来似乎就更自信了。因为之前曾经看过,“递归”是“动态规划”过程中使用的一种手段(当然,这句话也不准确。。。后面也会纠正它),这些特征表面,我使用的方法就算是“动态规划”了吧?
怀着将信将疑的态度,试着去查找了别人的文章,发现想了那么多层,最后自己还是没绕出“穷举”这个坑。/(ㄒoㄒ)/~~
4)分析纠正
前面说了几点“动态规划”特征的个人理解,很遗憾的是,那几点理解都不是严格地准确。这里我们来分析、纠正一下它们:
<1> 一个子问题如果属于一个大问题的一部分,就可以利用子问题的结果信息去计算大问题,避免重复计算
其实这并算是“动态规划”独有的特征。算法优化的本质就是为了“避免重复计算”。那么问题来了:都是在“避免重复计算”,不同的方法又有什么不同?它们避免重复计算的部分不一样吗?
而且,不是只有动态规划问题才能拆成子问题,例如贪心类型问题,把一段长的路拆成一步步的小路,把整体最优拆成计算局部最优,不也是可以大问题拆成小问题么?
<2>“递归”是“动态规划”过程中使用的一种手段
同理,并不是只有动态规划才会用到递归。“递归”说来,实际上是属于编程层面的东西,“函数自己调用自己”,是一种代码文字的书写、组织、表达方式,一种手段。而动态规划是一种思想,是可以脱离编程、程序而存在的,你可以在自己脑海中去运行这种思想。
2.动态规划
1)理论分析
正经的动态规划方法我参考学习了这篇文章:几个经典的动态规划问题
事实上,动态规划方法最重要的,是找到一个当前状态与前一状态对应的关系。因为在动态规划结构的问题中,当前状态会依赖于上一个状态(或认为是,当前的“可选”与之前的“选择”有关)。
如果拿数塔来说的话,就是我这一层能够选择哪些岔路,跟我之前走过、选择的路是有关系的。比如下图,我走了第一层的“9”,然后走了第二层的“15”,到了第三层我就只能走“6”、“8”这两种选择了,而不可能走“10”这个节点。这就说明我第三层的可选项跟我第一、二层的选择是有关系的。
如果是贪心结构的问题,比如下图,那当前的选择,就跟我之前的选择没有什么关系。比如我走到第二层,将要考虑如何选择第三层的路,不管我现在是处于H还是I,我的可选项都是一样(都是P、Q),我可以走到第三层的任意节点,因为H、I与第三层所有节点都是通的。每一步都随你瞎走,反正最后给你的选项都是一样的,至于路线是不是最优,那就是另说了。
如果把这种“可选项”理解成一种状态的话,那么动态规划结构的问题,每一层都是有很多种状态的,且不同的状态取决于上一层是的哪种状态。而贪心结构的问题,每一层的状态都只有一种。明白了动态规划问题存在着“状态”,并且状态与上个状态是存在着某种关系之后,我们去找出这种关系,将这种关系总结、归纳成“状态转换方程”,然后就是去编程实现了。动态规划问题的基本思路就是如此(当然这还是“正确的废话”,就算告诉你要去找状态转换方程,你也未必就能找出๑乛◡乛๑)
直接说吧,数塔问题的状态转换方程是这个:
dp(i, j) = Max( dp(i+1, j), dp(i+1, j+1 ) )+ f [ i ][ j ]
还需要把它翻译成人话:
为了便于说明,咱们先约定一下节点的表示方法吧:[ i, j ] ,它指的是第 i 层的第 j 个节点。比如下图[ 1, 1 ] 节点就是第1层、第1个节点,即“9”。
dp(i, j)函数的值,是[ i, j ] 这个节点,到最底层之后,所有可能路线的节点数值之和的最大值。还拿上图来说,比如dp(5, 1)其实就是[5, 1]这个节点(第5层、第1个)到最底层(第5层,也是它自己所在这一层)所有可能路线只有一条,即它自身“19”,所以dp[5, 1] = 19。又比如dp(4, 2),即[4, 2]这个节点,它到最底层(还是第5层)的所有路线有两条:7 -> 18、10 -> 18。数值和分别是25、28,28更大,所以最大值为28。因此dp(4, 2) = 28。
所有情况依次类推。
然后我们来解释一下这个状态转移方程式在说什么吧,dp(i, j) = Max( dp(i+1, j), dp(i+1, j+1 ) )+ f [ i ][ j ],意思是dp(i , j)这个的值,等于dp(i+1, j)、dp(i+1, j+1 ) 两者之间取大的,然后再加上 f [ i ][ j ]这个。
f [ i ][ j ]的含义是节点[i, j]的数值,比如上图[1, 1]节点,值是9,所以f[1][1 ]= 9。
再翻译得更通俗,就是你想要知道[i, j]这个节点到最底层的最大数值和,首先要知道[i+1, j]和[i+1, j+1]这两节点到最底层的最大数值和。知道之后,我们取大的,然后加上[i, j]节点的值,即f[i][j],就是我们要求的结果。
我们要求第1层、第1个节点到最底层的最大数值和,本质上就是要求dp(1, 1)。
2)实际操作
即使把公式翻译得通俗了,可能依然很绕口。我们就直接拿题目来说,一步步地来说明我们每一步怎么做。
想要知道根节点“9”到最底层的最大数值和,首先得知道第二层的“12”和“15”到最底层的最大数值和。我假设它们分别是M和N,假如“15”这条线路的数值和更大,即Max(M, N) = N,那么“9”到最底层的最大数值和,就是N + 9, 即下一层级的数值和 加上它自己。然后我们反推一步,按照第一步的类似原理,我们要知道“12”和“15”到最底层的最大数值和,就得知道第三层的“10”、”6“、”8“到最底层的最大数值和。
其中,连接“12”这个节点的是“10”、“6”。我假设这两条支线的数值和最大分别是O、P,假如”10“这条线路的更大,即Max(O、P) = O,那么“12”这条线到最底层的最大数值和就是 O + 12。
其中,连接“15”这个节点的是“6”、“8”。我假设这两条支线的数值和最大分别是P、Q(注意到,P就是上段话里的那个P,因为节点都是“6”,所以这个值是一样的),假如”6“这条线路的更大,即Max(P,Q) = P,那么“15”这条线到最底层的最大数值和就是 P + 12。
然后我们继续,往下一层推理、不断推理…
最后推理到了最底层,就没法继续往下推了,已经到了一个结束点。此时节点到最底层的数值和最大值,就是它自身。比如“19”,它到最底层线路的数值和最大就是它自己——19了。
现在我们明白了,要想知道第一层的信息,就必须知道第二层。 要想知道第二层的信息,又必须知道第三层…依次类推,只有最底层的信息,我们不需要由别的层来推导,而是可以直接得出。因此,最底层就是一个计算的起始点。
若从实现来说,确实也需要用到递归。最底层既是计算的起始点,同时又是递归的终止点、回溯点。
实现
说了那么多,我们就用程序来实现、验证吧。
1)数塔生成器
先来实现一个数塔生成器吧:
function createNumberTower (level) { let tower = new Array(); for (let i = 1; i < level + 1 ; i++) { let arr = new Array(i) for (let j = 0; j < arr.length; j++) { arr[j] = parseInt(Math.random() * 20) } tower.push(arr) } return tower; }
实际上就是生成一个二维数组,只是第二维的长度变化是有规律的,依次是1、2、3、4…这样自增。例如题目图里的数塔用数组表示就是[[9], [12, 15], [10, 6, 8], [2, 18, 9, 5], [19, 7, 10, 4, 16]]
参数是level,即要生成的数塔的层数。
这里我把每个节点设置为0~20的随机数,要想生成别的数,可以自行修改。
arr[j] = parseInt(Math.random() * 20)
2) 穷举法
虽然这个方法的效率低,但是既然已经思考出来了,那还是把它也一块实现了吧。穷举法的结果一般比较准确,正好可以用来验证我们的动态规划法的结果是否正确。
function findPath(tower) { let level = tower.length - 1; let prevPaths = tower.length === 2 ? [{path: [0], value: tower[0][0]}] : findPath(tower.slice(0, level)) let newPaths = new Array(); for (let path of prevPaths) { let lastNode = path.path[path.path.length - 1]; newPaths.push({ path: Array.from([...path.path,lastNode]), value: path.value + tower[level][lastNode] }) newPaths.push({ path: Array.from([...path.path,lastNode + 1]), value: path.value + tower[level][lastNode + 1] }) } return newPaths; }
结果将是一个数组,数组的每个元素是一条路径。里面包含两个字段信息:path、value。
path也是一个数组,用来表示路径每一层的节点索引。那题目图的数塔举例子,例如某条路径的path值是[0, 1, 2, 3, 4],那就说明第一层我走的是第1个元素(个数 = 索引 + 1嘛),第二层走的也是第2个元素,第三层走的是第3个元素,第四层走的是第4个元素,第五层走的是第5个元素。即9 -> 15 -> 8 -> 5 -> 16这一条路线。
value就是这一条路线所有节点的数值之和了。穷举了所用路线之后,按照value由大到小排序,就能知道哪一条路线数值之和最大了。
3)动态规划法
出乎意料,动态规划法的代码极其简单优雅,不过数行。
let fullTower = createNumberTower(30); // 生成数塔 function dp(i, j) { if (i === fullTower.length - 1) { // 起始点 return fullTower[i][j]; } else { return Math.max(dp(i+1, j), dp(i+1, j+1)) + fullTower[i][j]; } }
注意:为了编程方便,i,j 我用的是数组索引,即初始值是0。在这里dp(0 , 0)表示的才是对第1层第1个节点求解,而不是上文分析说明里的dp(1, 1)了。
4)验证
我们来验证一下,先把原题目求解一下。这里就不需要我们的数塔生成器了,直接把题目数塔写成一个二维数组
let fullTower = [[9], [12, 15], [10, 6, 8], [2, 18, 9, 5], [19, 7, 10, 4, 16]]
然后是穷举法:
findPath(fullTower).sort((prev, next) => next.value - prev.value) 输出: [{ path: [ 0, 0, 0, 1, 2 ], value: 59 }, { path: [ 0, 1, 1, 1, 2 ], value: 58 }, { path: [ 0, 0, 0, 1, 1 ], value: 56 }, { path: [ 0, 0, 1, 1, 2 ], value: 55 }, { path: [ 0, 1, 1, 1, 1 ], value: 55 }, { path: [ 0, 1, 2, 3, 4 ], value: 53 }, { path: [ 0, 0, 0, 0, 0 ], value: 52 }, { path: [ 0, 0, 1, 1, 1 ], value: 52 }, { path: [ 0, 1, 2, 2, 2 ], value: 51 }, { path: [ 0, 1, 1, 2, 2 ], value: 49 }, { path: [ 0, 0, 1, 2, 2 ], value: 46 }, { path: [ 0, 1, 2, 2, 3 ], value: 45 }, { path: [ 0, 1, 1, 2, 3 ], value: 43 }, { path: [ 0, 1, 2, 3, 3 ], value: 41 }, { path: [ 0, 0, 0, 0, 1 ], value: 40 }, { path: [ 0, 0, 1, 2, 3 ], value: 40 }]
可见,数值和最大的路径是[ 0, 0, 0, 1, 2 ]这种路线,也就是数塔图的9 - > 12 -> 10 -> 18 -> 10这种走法。该路线的和为59。
我们用动态规划来验证一下:
dp(0, 0) 输出:59
对比一下,发现与穷举法结果一致。证明我们的动态规划函数也是对的。为了简单性,这个dp函数里,我就不像穷举法一样,去列出具体的路径节点了。
我们再用数塔函数,随意生成一个数塔来试试。来生成一个层数更多的数塔吧,比如10层的:
let fullTower = createNumberTower(10); 输出:[ [ 16 ], [ 13, 0 ], [ 2, 14, 1 ], [ 10, 5, 17, 15 ], [ 14, 9, 10, 14, 16 ], [ 1, 7, 2, 18, 5, 19 ], [ 11, 0, 8, 5, 4, 11, 5 ], [ 13, 13, 11, 3, 2, 3, 9, 16 ], [ 3, 18, 3, 17, 13, 11, 18, 10, 18 ], [ 17, 2, 9, 7, 8, 2, 13, 8, 11, 5 ], [ 4, 5, 15, 0, 6, 9, 0, 12, 14, 19, 19 ], [ 7, 3, 17, 3, 3, 17, 13, 10, 11, 6, 5, 2 ], [ 14, 0, 5, 2, 9, 17, 10, 18, 6, 4, 6, 4, 19 ], [ 18, 9, 5, 12, 15, 8, 13, 0, 3, 2, 14, 17, 11, 15 ], [ 1, 6, 4, 18, 8, 3, 0, 1, 12, 12, 15, 6, 11, 7, 4 ], [ 7, 5, 4, 0, 1, 19, 10, 6, 18, 5, 16, 12, 2, 13, 10, 11 ], [ 1, 19, 11, 0, 15, 4, 15, 18, 0, 9, 14, 18, 4, 11, 1, 5, 7 ], [ 10, 8, 0, 8, 5, 11, 2, 1, 5, 16, 19, 16, 15, 14, 14, 4, 16, 15 ], [ 3, 19, 2, 4, 10, 2, 16, 1, 12, 10, 6, 5, 16, 15, 8, 7, 18, 10, 15 ], [ 2, 1, 10, 4, 12, 19, 11, 7, 16, 13, 0, 18, 14, 17, 0, 1, 13, 15, 5, 5 ], [ 16, 3, 1, 0, 1, 2, 5, 4, 8, 10, 17, 2, 18, 19, 6, 12, 8, 13, 10, 6, 16 ], [ 4, 15, 9, 6, 10, 17, 2, 11, 7, 17, 16, 2, 5, 2, 10, 15, 9, 15, 0, 10, 7, 8 ], [ 18, 18, 11, 11, 4, 1, 4, 13, 18, 13, 7, 18, 10, 15, 7, 5, 12, 18, 10, 15, 12, 0, 9 ], [ 8, 8, 11, 18, 6, 4, 17, 15, 7, 2, 7, 12, 12, 11, 9, 7, 17, 15, 5, 3, 12, 10, 5, 7 ], [ 14, 2, 5, 16, 1, 6, 7, 6, 10, 15, 6, 1, 6, 3, 16, 4, 8, 0, 16, 5, 11, 10, 9, 16, 1 ], [ 13, 0, 13, 3, 19, 11, 3, 9, 5, 5, 11, 11, 15, 7, 9, 3, 11, 2, 0, 2, 11, 19, 5, 1, 18, 12 ], [ 10, 10, 16, 12, 14, 3, 19, 11, 9, 0, 10, 6, 7, 0, 6, 14, 5, 1, 11, 10, 6, 13, 5, 5, 3, 14, 5 ], [ 3, 5, 4, 14, 13, 15, 1, 8, 18, 14, 5, 16, 11, 5, 10, 6, 15, 10, 10, 6 ], [ 10, 18, 7, 17, 1, 10, 17, 9, 3, 17, 16, 10, 3, 1, 17, 1, 12, 6, 9, 19, 13, 18, 14, 0, 5, 5, 0, 19, 17 ], [ 0, 17, 18, 16, 2, 5, 9, 13, 15, 0, 18, 13, 8, 0, 0, 15, 16, 4, 12, 10, 12, 8, 15, 16, 2, 9, 1, 4, 7, 15 ] ]
穷举法的结果是:
findPath(fullTower).sort((prev, next) => next.value - prev.value) 输出:[ { path: [ 0, 1, 2, 3, 4, 4, 5, 6, 7, 8 ], value: 126 }, { path: [ 0, 1, 1, 1, 2, 3, 3, 4, 5, 5 ], value: 123 }, { path: [ 0, 1, 1, 2, 3, 4, 5, 6, 7, 8 ], value: 122 } ....]
动态规划法给出的结果是:
dp(0, 0) 输出:126
可见也是一致的。多次尝试,结果都是一致,这样就验证了我们想法的正确性。
-
OR Paper Weekly(一) | 用机器学习生成列生成的列,元启发式算法=动物世界?看OR68年发文数据,哪国位居...
2022-02-26 15:43:35辅之以科普/点评/吐槽的方式,让大家随时了解最新的科研动态。欢迎大家一起来 欣赏优质文章,学习脑洞文章,鄙视灌水文章。 精选论文 (一) 论文题目:Machine-Learning–Based Column Selection for C作者:王源,徐思坤,陈贤邦
OR Paper Weekly 栏目将会从运筹学顶级期刊上选择一部分有趣的文章,对这些文章的主要研究内容进行一个概述/点评。OR Paper Weekly 的特点是 不做大而全的照搬,也未必都只选择优质的文章,而是精选一部分有趣的文章。辅之以科普/点评/吐槽的方式,让大家随时了解最新的科研动态。欢迎大家一起来 欣赏优质文章,学习脑洞文章,鄙视灌水文章。
精选论文 (一)
论文题目:Machine-Learning–Based Column Selection for Column Generation
期刊:Transportation Science
发表年份:2020年
作者:Mouad Morabit, Guy Desaulniers, Andrea Lodia
原文链接:
https://pubsonline.informs.org/doi/abs/10.1287/trsc.2021.1045
摘要:列生成算法是常见的求解大规模整数规划问题的方法。本文提出了一种基于机器学习来加速列生成求解的新方法。该方法的特点是 在每一步列生成算法的迭代中,通过机器学习训练出一个模型来生成新的列(选择一个变量的子集添加到主问题构成新的一列)。该方法通过选择最合适的列来降低主问题的计算时间。我们在带有时间窗约束的车辆路径规划问题和组环问题上验证了该方法。实验结果表明该方法可以降低30%的计算时间。
文章亮点/点评:传统的列生成算法是通过极大化 reduced cost 的子问题来生成新的列的。从全局的视角来看 这样添加列的方式会比较的greedy,因此传统列生成算法的效率会受到较大影响。本文通过机器学习来生成新的列,将机器学习和传统的运筹优化算法相结合,以基于数据的机器学习模型来解决传统运筹学中一些难以解决的问题已经成为了目前的研究趋势。
精选论文 (二)
论文题目:Metaheuristics “In the Large”
期刊:European Journal of Operational Research
发表年份:2021年
作者:Jerry Swan, Steven Adriaensen, Alexander E. I. Brownlee c, Kevin Hammond d, Colin G. Johnson, et al.
原文链接:
https://www.sciencedirect.com/science/article/pii/S0377221721004707
摘要:在过去的几十年里涌现出大量关于元启发式算法(Metaheuristics)的各类改进,可以说元启发式算法早已成为了优化领域的显学。为了避免元启发式算法缺乏系统性和缺乏可重复性的问题,我们迫切需要强大的计算设备来支持新方法的开发,分析和比较。为此,我们提出了 Metaheuristics “In the Large” 项目,以实现我们的愿景。该项目的愿景是:实现真正可扩展的算法模版,无需对算法进行修改即可重用。为特定领域的知识(domain knowledge)的融入提供通用的白盒问题描述,以及可以远程访问的框架,组建和问题。相信这一举措会提高算法的可重复性,加速该领域的进一步发展。我们认为通过这种基础施设的建设,可以促使该领域朝着更水平的科学研究加速前进。在本文中我们描述了我们的期望,并且展示了所有元启发采用的通用协议,以帮助进一步释放该领域的潜力,简化对元启发式算法设计空间的探索。
文章亮点/点评:诚如本篇论文所述,元启发式算法在过去的几十年中出(guan)现(shui)了大量的论文。各种各样的生物在元启发式算法中粉墨登场,例如蛙跳算法、猫群算法、蟑螂算法,天牛算法,帝王蝶算法,头脑风暴算法,黏菌觅食算法,植物胞群算法,树种优化算法。由此可见,动物都不够用的了,元启发式算法已经将触手伸到了昆虫,植物和微生物。元启发式算法有着非常好的应用场景,例如在多目标优化,数据驱动优化和黑箱优化等场景下可以和数学优化算法形成互补。希望通过Metaheuristics “In the Large” 项目能促使元启发式算法的研究进一步的向着健康的方向发展。
精选论文 (三)
论文题目:Operations Research: Topics, Impact, and Trends from 1952–2019
期刊:Operations Research
发表年份:2021年
作者:Mouad Morabit, Guy Desaulniers, Andrea Lodia
原文链接:
https://pubsonline.informs.org/doi/abs/10.1287/opre.2021.2139
摘要:本文回顾了 Operations Research 期刊近68年来的出版成果,揭示了其文章、作者的变化及其随时间的影响,以及这些变化对研究人员和从业者的影响。Operations Research 自1952年创立到2019年,总共接收了5440篇文章。
通过对论文作者的数据进行分析可得以下结论:1美国、加拿大和英国分列发文的前三位,其中美国以4058篇名列第一,并且美国的发文量超过了第2名到第10名的总和。中国(大陆)以125篇名列第7名。2 库存问题是北美、亚洲和中东国家最热门的研究问题,而欧洲国家则专注于调度问题。
通过对研究内容的数据进行分析可得以下结论:1 数学规划是最常见的研究方法,2 库存问题是研究最多的问题,3近年来与定价有关的研究正在显著增加。
通过对研究趋势的数据进行分析可得以下结论:1 近十年来动态规划是最常用的方法,2 近十年来定价是研究最多的问题。由此可知动态规划和定价问题在未来依然会是研究热点。
文章亮点/点评:
关于本文的详细解读 可以参考之前我们的一篇文章:【主编推荐】哪个国家在 Operations Research 上发文最多?(数据解析OR创刊68年来发文趋势)
精选论文 (四)
论文题目:From Predictive to Prescriptive Analytics
期刊:Management Science
发表年份:2018年
作者:Dimitris Bertsimas, Nathan Kallusb
原文链接:
https://pubsonline.informs.org/doi/abs/10.1287/mnsc.2018.3253
摘要:我们构造了一个结合了机器学习、运筹与管理科学以及一些特殊的方法的框架,通过这个框架,我们可以利用数据为运筹与管理科学中的问题推荐(Prescribe)最优决策。与其他数据驱动优化的研究不同的是,我们考虑的数据不仅包含了对成本和收益有直接影响的观测量,还包含了问题相关的主要辅助变量的观测值。本文考虑的主要问题是一个基于不完美样本的条件随机优化问题,其中问题相关的联合概率分布未知。我们展示了我们的方法可广泛应用于非常多的决策问题,证明了这些方法在计算上易解并且在一般条件下渐进最优,这些结果即使在不满足独立同分布的数据以及删失数据下仍然成立。我们将这套方法拓展至某些决策变量(例如价格)与不确定性有未知因果关系的场景。我们定义了指导性系数(Coefficient of Prescriptiveness)P,用以从运营的角度度量数据的指导性内容以及策略的有效性。我们以一个年均发货量超过十亿单的大型媒体公司分销部门所面对的库存管理问题为例来展示这套方法。我们利用内部数据以及从IMDb、烂番茄和谷歌上收集得到的公共数据对他们的运营决策所作出的指导性建议提升了他们运营的基准表现。具体而言,利用我们的方法,我们收集到的数据在指导性系数的度量下为他们的运营能力带来了88%的提升。
文章亮点/点评:从描述性分析(Descriptive analysis)中的收集、整理、统计,到预测性分析(Predictive analysis)中更为复杂的挖掘、聚类、预测,我们对数据的了解越来越深。现在我们不仅需要理解数据,更要利用数据来做出更智能的决策,这就是本文所描述的指导性分析(Prescriptive Analysis)。数据驱动决策是目前绝大多数商业模式迈向智能与高效的关键,而实现这样的跨越需要背后完备的概率统计和优化理论的支撑。本文从概述利用数据推荐最优决策的框架出发,通过研究一类特殊的随机优化问题,以及它在经典的库存管理问题中的应用,展示了数据驱动决策在运营效率提升上的强大之处。虽然 Prescriptive Analysis一词有新瓶装旧酒之嫌,但无论是为了了解数据驱动决策的框架还是为了学习它在具体运营问题中的应用,本文都是一篇值得一读的经典。
精选论文 (五)
论文题目:Actor-Critic–Like Stochastic Adaptive Search for Continuous Simulation Optimization
期刊:Operation Research
发表年份:2021年
作者:Qi Zhang, Jiaqiao Hua
原文链接:
https://pubsonline.informs.org/doi/abs/10.1287/opre.2021.2214
摘要:我们为求解一类带Lipschitz连续性质的仿真优化问题提出了一个随机搜索算法。该算法从一个参数化的解空间概率分布中抽取候选解,并通过一个基于所谓的收敛球方法的异步学习过程衡量抽取的样本点(候选解)的表现。这个算法最大的特点在于它在完整地保留了此前仿真的信息的同时采用了一个近似架构来利用对目标函数的知识来提升解的质量。算法的每一步同时调整参数化的解空间概率分布以及目标函数的近似估计,这与强化学习中的Actor-Critic架构非常类似。我们提出了算法表现的有限时间概率边界,并展示了在每次迭代中仅有一个仿真观测值的情况下的全局收敛性。我们的实证研究表明该算法很有前景,并有可能在效率以及可靠性上超越目前的一些搜索方法。
文章亮点/点评:本文所提出的随机搜索算法,在传统随机优化的思路,即近似目标函数的基础上,赋予解空间一个先验并用随机抽样的方法去学习解空间的分布,来达到提升算法表现的目的。文章算法的思路非常清晰,并且有一定的启发性。与强化学习中Actor-Critic算法的类比非常抓人眼球。这也从侧面说明了一篇顶刊文章不仅需要方法论上的创新,还需要讲好一个故事,有一个抓人眼球的点。
精选论文 (六)
论文题目:Feature-driven Economic Improvement for Network-constrained Unit Commitment: A Closed-loop Predict-and-optimize Framework
期刊:IEEE Transactions on Power Systems
发表年份:2021年
作者:Xianbang Chen, Yafei Yang, Yikui Liu, Lei Wu
原文链接:
https://www.researchgate.net/profile/Xianbang-Chen/publication/356189174_Feature-driven_Economic_Improvement_for_Network-constrained_Unit_Commitment_A_Closed-loop_Predict-and-optimize_Framework/links/61906a2961f098772093c4d3/Feature-driven-Economic-Improvement-for-Network-constrained-Unit-Commitment-A-Closed-loop-Predict-and-optimize-Framework.pdf
Abstract:
机组组合(unit commitment, UC)是电力市场结算的重要环节。一般来说,UC会在开环型预测-再优化(Open-loop predict-then-optimize, O-PO)框架下执行:先以统计精度为导向预测可再生能源出力,再以系统运行成本最小为目标优化UC,预测与优化环节呈开环状态。但由于UC模型含有大量非线性约束,因此统计上更精确的预测值未必会诱导出经济性更好的UC决策。针对该问题,本文提出了一种闭环型预测-优化(Closed-loop predict-and-optimize, C-PO)框架以提高UC经济性。C-PO以UC经济性为导向进行预测,同时考虑了UC模型的约束与优化目标,使得预测与优化环节呈闭环状态。为提高C-PO的实用性,设计了拉格朗日分解算法,以加速C-PO预测器的训练过程。最后,基于IEEE RTS 24-Bus系统和5655-Bus系统的算例得到两个结论:C-PO框架能够有效提高UC经济性,并具有较好的实用性。
文章亮点/点评:
1. 生成以UC经济性为导向的可再生能源预测值,C-PO能够有效提高UC经济性。
2. 无需替换已有的预测工具,C-PO与当前普遍使用的O-PO框架兼容。
3. 拉格朗日分解算法能有效加速预测器的训练,C-PO具备应用于大型系统的潜力。
-
动态规划常见类型总结
2019-03-26 23:55:28严格来说,递推不属于动态规划问题,因为动态规划不仅有递推过程,还要有决策(即取最优),但广义的动态规划是可以包含递推的,递推是一类简单的、特殊的动态规划,毕竟动态规划与递推密不可分。动态规划类型主要... -
动态规划算法经典案例
2019-04-29 19:20:10动态规划算法是从暴力搜索算法优化过来的,如果我们不清楚暴力搜索的过程,就难以理解动态规划的实现,当我们了解了动态规划算法的基本原理的文字概述,实现条件之后,这时可能并不是太理解这种思想,去面对实际问题... -
c++动态规划经典算例
2019-06-05 09:46:07动态规划算法通常用于求解具有某种最优性质的问题。在这类问题中,可能会有许多可行解。每一个解都对应于一个值,我们希望找到具有最优值的解。动态规划算法与分治法类似,其基本思想也是将待求解问题分解成若干... -
实例详解贪心、分治、回溯以及与动态规划的比较分析
2019-07-31 11:35:19上一章已经详细地分析过动态规划,接下来还是结合具体问题,感受一下这些算法是怎么工作的,是如何解决问题的,再问体体会这些算法的本质。我觉得比单纯记忆原理和定义更有价值。 贪心算法 1. 如何理解贪心算法? ... -
动态规划C++
2018-11-09 10:56:46引入题目 ...分析:可以使用,暴力搜索方法、记忆搜索方法、动态规划方法、状态继续简化后的动态规划方法。 暴力搜索方法 arr = {5,10,25,1},aim=1000. 1.用0张5元货币,让[10,25,1]组成剩下的1... -
动态规划
2016-10-28 17:01:50动态规划(DP)[1]通过分解成子问题解决了给定复杂的问题,并存储子问题的结果,以避免再次计算相同的结果。我们通过下面这个问题来说明这两个重要属性:重叠子问题和最优子结构。 -
立体匹配---动态规划
2017-05-22 10:58:18近来研究立体匹配,从入门开始,先学习一些基本的算法思想。 立体匹配算法中,全局匹配是一个很重要的部分,利用...全局匹配算法通过构建全局能量函数,然后通过优化方法最小化全局能量函数以求得致密视差图。 -
动态规划详解(第三讲)——矩阵链相乘
2018-07-30 21:55:07的列数一定等于矩阵 M i + 1 M i + 1 M_{i+1} 的行数( 1 ≤ i < n 1 ≤ i < n 1\le i ),这是由矩阵乘法的定义决定的。 因此,对于一个矩阵链,我们指定每个矩阵的行数和最右面矩阵 M n M n M_n 的列数就可以... -
动态规划从理论到实践-深入理解贪心/分治/回溯/动态规划的算法思想
2019-04-27 22:49:06动态规划:一个模型(多阶段决策最优解模型),三个特征(最优子结构、无后效性和重复子问题),适用于求解最优问题。本文主要参考资料是王争的《数据结构与算法之美》,数据结构和算法作为程序员的... -
巧解动态规划问题
2020-02-05 09:11:40本文是对动态规划问题的总结,看完之后会对类似的题目有所了解。 -
动态规划算法题总结
2018-08-09 17:59:20动态规划的基本思想 动态规划(Dynamic Programming,简称DP),虽然抽象后进行求解的思路并不复杂,但具体的形式千差万别,找出问题的子结构以及通过子结构重新构造最优解的过程很难统一。 在做这些题之前,... -
分治、贪心、回溯、动态规划四大算法
2020-02-14 16:27:00文章目录 本文转于这,谢谢原作者 前言 根据自己对贪心算法、分治算法、回溯算法、动态规划四种算法思想的理解对其分别做一个引入和介绍。 参考极客时间王争老师的数据结构与算... -
算法导论读书笔记(15)动态规划
2017-06-16 17:40:17其速度比动态规划方法快得多,但是,我们并不总能简单地判断出贪心算法是否有效 3. 摊还分析并不是通过分别分析每个操作的实际代价的界来分析操作序列的代价的界,而是直接分析序列整体的实际代价的界,好处是,... -
【planning in frenet Frame】局部路径、速度解耦规划的一般方法
2022-04-15 08:59:22系列文章目录 提示:这里可以添加系列文章的所有文章的目录,目录需要自己手动添加 ...【速度规划】在frenet坐标系的ST图上进行最优路径path速度规划,生成平滑速度轨迹(1)连续空间的离散化,设置栅格 -
C++ 数字三角形(动态规划)
2019-01-25 10:29:08一、题目 有三种方法实现: 这里就求权值之和...动态规划 思路和自底向上、自顶向下一样,只不过开多一个数组opt存储计算过的到达a[i][j]的权值之和最大的数 二、思路: 方法1:自底向上 1.从最后一行出... -
【动态规划 & 备忘录】矩阵连乘问题(C++)
2020-05-28 15:19:22自顶向下 递归 的动态规划2. 自底向上 非递归 的动态规划3. 运行结果展示 一、关于矩阵连乘 1. 问题描述: 给定n个矩阵:A1,A2,…,An,其中Ai与Ai+1是可乘的(i=1,2…,n-1)。确定计算矩阵连乘积的计算次序,使得... -
动态规划C++实现--最小编辑代价
2018-05-15 20:41:18如果 str1 的长度为 M,str2的长度为 N,采用动态规划的方法,时间复杂度为 O(M*N) ,额外的空间复杂度为O(M*N) 首先生成动态矩阵 dp[M+1][N+1],其中 dp[i][j] 表示 st1[0,...,i-1] 编辑成 str2[0,...,j-1]的... -
动态规划100例
2015-05-11 11:57:14更多"动态规划100例"相关资料请点击这里 一 资源问题.....................................................................................................................................4 资源问 -
五类常见算法小记 (递归与分治,动态规划,贪心,回溯,分支界限法)
2014-05-31 23:15:275类算法小结: 递归与分治法, 动态规划, 贪心算法, 回溯法, 分支界限法 -
立体匹配中的全局匹配(一)动态规划笔记
2016-09-19 16:17:38近来研究立体匹配,从入门开始,先学习一些基本的算法思想...全局匹配算法一般有动态规划、置信传播、模拟退火、图割法、遗传学等,这里首先介绍动态规划,也是从一些论文中提取的思想,可能有不对的地方,望指正。动态 -
算法基础-->贪心和动态规划
2017-09-05 23:23:56本篇博文将详细总结贪心和动态规划部分,贪心和动态规划是非常难以理解和掌握的,但是在笔试面试中经常遇到,关键还是要理解和掌握其思想,然后就是多刷刷相关一些算法题就不难了。本篇将会大篇幅总结其算法思想。... -
【算法设计与分析】矩阵连乘问题(动态规划)
2019-04-08 18:02:24一、背景介绍 ...原因:矩阵乘法满足结合律; 这种计算次序可以用加括号的方式来确定。 完全加括号 若一个矩阵连乘积的计算次序完全确定,也就是说该连乘积已完全加括号; 则可以依此次序反复调用2个矩阵...