利用MATLAB求解离散动态规划问题.pdf

利用MATLAB求解离散动态规划问题.pdf

《利用MATLAB求解离散动态规划问题.pdf》由会员分享，可在线阅读，更多相关《利用MATLAB求解离散动态规划问题.pdf（28页珍藏版）》请在冰豆网上搜索。

利用MATLAB求解离散动态规划问题利用MATLAB求解离散动态规划问题陈俊先中山大学岭南学院2016年8月30日1/28利用MATLAB求解离散动态规划问题一、简介目前MATLAB内置的最优化工具包不能很好的处理多期最优化问题，在求解动态规划问题上并不非常适用目标：

利用MATLAB实现离散动态规划问题的求解，寻求并设计编写一种通用的求解算法和相应的函数试图利用MATLAB并用最简单易懂的方法解决讲义中所有例题与练习2/28利用MATLAB求解离散动态规划问题二、穷举法：

算例1例2.1动态规划讲义第35页习题2，设国家拨给60万元投资，供四个工程扩建使用，每个工程扩建后的利润与投资额的大小有关，投资后的利润函数如下表所示。

问如何投资使得总利润最大化？

利润投资0102030405060v1（x）0205065808585v2（x）0204050556065v3（x）0256085100110115v4（x）0254050606570穷举所有可能的结果，找出其中最优值及对应的策略3/28利用MATLAB求解离散动态规划问题二、穷举法：

算例1-算法与程序对于上述问题，在穷举求解时应注意投资额不多于60万元的约束。

稍加分析便知道60万元必须全部用完没有剩余，否则结果一定不会是最优解。

因此可以写出以下程序:

输入利润函数：

V1=0,20,50,65,80,85,85;V2=0,20,40,50,55,60,65;V3=0,25,60,85,100,110,115;V4=0,25,40,50,60,65,70;4/28利用MATLAB求解离散动态规划问题二、穷举法：

算例1-算法与程序（续）m=1;J=zeros（84,1）;U=zeros（84,4）;fori=1:

7forj=1:

7fork=1:

7forl=1:

7ifi+k+j+l=10J（m）=V1（i）+V2（j）+V3（k）+V4（l）;U（m,1）=i;U（m,2）=j;U（m,3）=k;U（m,4）=l;m=m+1;elsecontinue;end.5/28利用MATLAB求解离散动态规划问题二、穷举法：

算例1-运行结果得到最优值与相应的最优策略：

OptV=max（J（1:

end）OptU=（U（find（J=OptV）,:

）1）10运行结果：

OptV=160OptU=20,0,30,10由于计算量很小，整个过程的运算速度非常快，从MATLAB自带的计时器来看只用了大约0.009秒。

6/28利用MATLAB求解离散动态规划问题二、穷举法：

算例2例6.1最短路问题。

如下图所示，其中小圆圈称为点，两点间的连线称为弧，弧上的数字称为弧长，箭头表示前进方向。

要寻找一条从起点A到终点E的连通弧，使其总弧长最短。

7/28利用MATLAB求解离散动态规划问题二、穷举法：

算例2在开始编程求解之前需要先简单分析一下题目的情景，以便简化求解流程。

稍加分析就能发现，其实这个问题只需要做三次决策，分别在A、B和C阶段，行进到D阶段以后下一步自然就只能到达E点，也就是终点。

根据以上分析以及最短路问题的特点，将各结点用数字标明如下图：

8/28利用MATLAB求解离散动态规划问题二、穷举法：

算例2首先根据最短路问题中状态转移的特点，写出带权的状态转移矩阵如下：

1234567891011153213638764685356337848393102119/28利用MATLAB求解离散动态规划问题二、穷举法：

算例2-算法与程序在上面的带权的状态转移矩阵中，若两个结点间不能相互到达或按着题给示意图的方向到达则留空，若能够到达，则表格中的数值代表弧长。

注意这里的到达关系是单项的，只能顺着题给示意图的方向行进。

对应上面的状态转移矩阵，在MATLAB中对应初始化，其中inf代表无穷大数，表示两个结点之间不能相互到达或不能按题给示意图方向到达：

inf,5,3,inf,inf,inf,inf,inf,inf,inf,infinf,inf,inf,1,3,6,inf,inf,inf,inf,infinf,inf,inf,inf,8,7,6,inf,inf,inf,infinf,inf,inf,inf,inf,inf,inf,6,8,inf,infinf,inf,inf,inf,inf,inf,inf,3,5,inf,infR=inf,inf,inf,inf,inf,inf,inf,inf,3,3,infinf,inf,inf,inf,inf,inf,inf,inf,8,4,infinf,inf,inf,inf,inf,inf,inf,inf,inf,inf,3inf,inf,inf,inf,inf,inf,inf,inf,inf,inf,2inf,inf,inf,inf,inf,inf,inf,inf,inf,inf,2inf,inf,inf,inf,inf,inf,inf,inf,inf,inf,inf;10/28利用MATLAB求解离散动态规划问题二、穷举法：

算例2-算法与程序（续）n=size（R,1）;Path=zeros（50,5）;Path（1:

50,1）=1;Path（1:

50,5）=n;D=infones（50,1）;m=1;fori=1:

n;ifR（1,i）=infforj=1:

nifR（i,j）=inffork=1:

nifR（j,k）=infPath（m,2:

4）=i,j,k;D（m）=0;D（m）=D（m）+R（1,i）+R（i,j）+R（j,k）+R（k,n）;m=m+1;end.11/28利用MATLAB求解离散动态规划问题二、穷举法：

算例2-运行结果最后照例在穷举结果中找最优值：

OptD=min（D（1:

end）OptPath=Path（find（D=OptD）,:

）运行结果：

OptD=14OptPath=1,2,5,8,11翻译过来就是AB1C2D1E，这与我们传统求解方法得到的结果相一致。

整个过程依旧非常迅速，耗时约0.008秒。

最后还顺带得到所有可能的路径为12条。

12/28利用MATLAB求解离散动态规划问题二、穷举法：

优缺点优点：

程序简单，容易理解；在阶段数量与状态数量并不很大时运算速度快缺点：

并未运用动态规划思想，对不同问题通用性低，程序需要较大改动13/28利用MATLAB求解离散动态规划问题三、动态规划求解函数functionpopt,fval=dynprog（x,DecisFun,ObjFun,TransFun）输入参数：

x是状态变量，一列代表一个阶段的状态函数DecisFun（k,x,u）描述可行决策集ObjFun（k,s,u）描述阶段值TransFun（k,s,u）描述状态转移方程输出变量：

popt=序号组，状态变量组，最优决策组，阶段值组fval也是一个列向量，表示最优值14/28利用MATLAB求解离散动态规划问题三、动态规划求解函数-算例1求解过程分析1.阶段变量k：

将问题按工程分为4个阶段，即k=1,2,3,42.状态变量sk：

状态变量sk表示第k个工程至第4个工程可用的投资总额，显然有s1=6015/28利用MATLAB求解离散动态规划问题三、动态规划求解函数-算例1求解过程分析（续）3.决策变量uk：

允许决策集：

0uksk，特别地，当k=4，uk=sk因此可写出以下可行决策集函数：

functionu=DecisFun（k,s,u）ifk=4u=s;elseu=0:

10:

s;end16/28利用MATLAB求解离散动态规划问题三、动态规划求解函数-算例1求解过程分析（续）4.状态转移方程Tk：

显然有：

sk+1=sku因此可写出以下状态转移函数：

functionsnext=TransFun（k,s,u）snext=su;17/28利用MATLAB求解离散动态规划问题三、动态规划求解函数-算例1求解过程分析（续）5.阶段值Vk：

第k阶段的阶段值表示投资给第k个工程uk时的利润因此可写出以下阶段值函数：

functionV=ObjFun（k,s,u）w=0000202025255040604065508550805510060856011065856511570w=w;V=w（u/10+1,k）;18/28利用MATLAB求解离散动态规划问题三、动态规划求解函数-算例1求解过程分析（续）6.各阶段状态变量可能取值：

由已知：

s1=60s2=0,10,20,30,40,50,60s3=0,10,20,30,40,50,60s4=0,10,20,30,40,50,60因此状态变量可能的取值可以描述为：

s=nanones（7,4）;s（1,1）=60;s（:

2）=01020304050600;s（:

3）=01020304050600;s（:

4）=01020304050600;19/28利用MATLAB求解离散动态规划问题三、动态规划求解函数-算例1求解过程分析（续）调用动态规划求解函数：

popt,fval=dynprog（s,0DecisFun0,0ObjFun0,0TransFun0）popt=16020502400034030854101025fval=160对于fval，因为之前把阶段值转变为负值，这个结果是-160，而对于我们的问题中应是160。

运用动态规划求解函数求解总共耗时0.017秒，相比穷举法耗时更多，因为在阶段数和状态数较少时，求解函数不能体现出优势，调用函数的过程涉及各种参数的传递和运算，这是相对比较耗时的。

20/28利用MATLAB求解离散动态规划问题三、动态规划求解函数-算例3例18.1机器负荷分配问题。

某机器可以在高低两种负荷下生产。

如果年初投入的完好机器数为u1，那么在高负荷生产时，年终的完好机器数为0.7u1（系数0.7称为机器的完好率），年产量为8u1；在低负荷生产时，完好率为0.9，年产量为5u1。

假设开始生产时期有1000台完好机器。

要指定一个五年计划，在每年年初决定重新分配在两种负荷下生产的完好机器数量，使五年内产品的总产量最大。

运用动态规划求解函数求解这个问题只需要按照上面的形式套用即可，分析如下：

21/28利用MATLAB求解离散动态规划问题三、动态规划求解函数-算例3求解过程分析1.阶段变量k：

将问题分成五个阶段，即k=1,2,3,4,52.状态变量sk：

状态变量sk表示第k个工程至第4个工程可用的投资总额，显然有s1=100022/28利用MATLAB求解离散动态规划问题三、动态规划求解函数-算例3求解过程分析（续）3.决策变量uk：

允许决策集：

0uksk因此可写出以下可行决策集函数：

functionu=DecisFun（k,s,u）u=0:

s;这里以1为间隔，当然可以分的更细，但运行时间较长23/28利用MATLAB求解离散动态规划问题三、动态规划求解函数-算例3求解过程分析（续）4.状态转移方程Tk：

显然有：

sk+1=sku因此可写出以下状态转移函数：

functionsnext=TransFun（k,s,u）snext=0.9s0.2u;24/28利用MATLAB求解离散动态规划问题三、动态规划求解函数-算例3求解过程分析（续）5.阶段值Vk：

functionV=ObjFun（k,s,u）V=5s3u;只要注意这里是负值即可，若是求最小，则还是正25/28利用MATLAB求解离散动态规划问题三、动态规划求解函数-算例3求解过程分析（续）6.各阶段状态变量可能取值：

由已知：

s1=1000s2=0,1,2,.,998,999,1000s3=0,1,2,.,998,999,1000s4=0,1,2,.,998,999,1000s5=0,1,2,.,998,999,1000因此状态变量可能的取值可以描述为：

s=nanones（1001,5）;s（1,1）=10000;s（:

2）=0:

10000;s（:

3）=0:

10000;s（:

4）=0:

10000;s（:

5）=0:

10000;26/28利用MATLAB求解离