动态规划编程
- 格式：doc
- 大小：93.00 KB
- 文档页数：9
动态规划编程动态规划是解决多阶段决策过程最优化问题的一种方法．该方法是由美国数学家贝尔曼(R．Bellman)等人在2O世纪50年代初提出的．他们针对多阶段决策问题的特点，提出了解决这类问题的最优化原理，并成功地解决了生产管理、资源分配等方面的许多实际题，从而建立了运筹学的一个新分支——动态规划．动态规划是现代企业管理中的一种重要决策方法，可用于解决最优路径、资源分配、生产计划与库存、投资、装载、排序等问题，还可用于生产过程的最优控制等．由于它有独特的解题思路，因而在处理某些优化问题时，比线性规划或非线性规划方法更有效．而Matlab是一个功能强大的用于基于矩阵运算的强大数值计算软件，将Matlab语言应用到动态规划中去，对实际问题进行程序设计和计算，可以达到计算简便的目的．一、动态规划基本概念使用动态规划方法解决多阶段决策问题，首先要将实际问题写成动态规划模型，此时要用到以下概念：1)阶段将所给问题的过程，按时间或空间特征分解成若干互相联系的阶段，以便按次序去求解每阶段的解，每个阶段就是一个子问题，常用字母k表示阶段变量．2)状态各阶段开始时的客观条件叫做状态．描述各阶段状态的变量称为状态变量，常用sk表示第k阶段的状态变量．状态变量sk 的取值集合称为状态集合，用sk表示．3)决策当各段的状态取定以后，就可以作出不同的决策(或选择)，从而确定下一阶段的状态，这种决定称为决策．表示决策的变量称为决策变量，常用uk (sk)表示第k阶段当状态为sk时的决策变量．在实际问题中，决策变量的取值往往限制在一定范围内，称此范围为允许决策集合，常用Dk (sk)表示第k阶段从状态sk出发的允许决策集合，显然有“uk∈Dk(sk)．4)策略一个由每个阶段的决策按顺序组成的集合称为策略，用p表示，即p(s1)={u1(s1),u2(s2),.......,un(sn),}。
一个n阶段决策过程，从1到n叫作问题的原过程．对于任意给定的k(1≤k≤n)，从第k阶段状态sk 到第n阶段状态sn的过程称为原过程的一个后部子过程．后部子过程的策略记为pk (sk)={uk(sk),uk+1(sk+1),......,un(sn)}，在实际问题中，可供选择的策略有一定的范围，此范围成为允许策略集合。
允许策略集合中达到最优效果的策略成为最优策略5)状态转移动态规划中本阶段往往是上一阶段状态和上一阶段的决策进行综合的结果．如果给定了第k段的状态sk ，且该阶段决策为uk(sk)，则第k+1段的状态sk+1也就完全确定．它们的关系可表示为：s k+1=Tk(sk,uk)由于上式表示了由k阶段到k+1阶段的状态转移规律，所以称该式为状态转移方程．6)指标函数用于衡量所选定策略优劣的数量指标称为指标函数．一个n阶段决策过程，从1到n叫作问题的原过程．对于任意一个给定的k(1≤k≤n)，从第k阶段到第n阶段的过程称为原过程的一个后部子过程。
V1,n (s1,p1,n)表示初始状态为s1采用策略p1,n时原过程的指标函数值。
而V k,n (sk,pk,n)表示在第k阶段，状态为sk采用策略pk,n时后部子过程的指标函数值．最优指标函数记为fk (sk)，它表示从第k阶段状态sk采用最优策略pk,n到过程终止时的最佳效益值．fk (sk)与Vk,n(sk,pk,n)间的关系为：f k (sk)=Vk,n(sk,pk,n)=optimize Vk,n(sk,pk,n)当k=1时，f1(s1)就是从初始状态s1到全过程结束的整体最优函数．二、动态规划基本思路1)将多阶段决策过程划分阶段，恰当地选择状态变量、决策变量以定义最优指标函数，从而把问题化成一族同类型的子问题，然后逐个求解．2)求解时从边界条件开始，逆序过程行进，逐段递推寻优．在每一个子问题求解时，都要使用它前面已求出的子问题的最优结果．最后一个子问题的最优解，就是整个问题的最优解．3)动态规划方法是既将当前一段与未来各段分开，又把当前效益和未来效益结合起来考虑的一种最优化方法，因此每段的最优决策选取是从全局考虑的，与该段的最优选择一般是不同的．三、动态规划函数使用说明由于我们的目的是使用动态规划解题，而不是要我们直接编写动态规划的MATLAB程序，那我们下面直接给出一个现成的动态规划的MATLAB源代码，这里不讨论它到底是是如何运行的，我们只是说明下该函数如何使用1.function [p_opt,fval]=dynprog(x,DecisFun,ObjFun,TransFun)2.% 自由始端和终端的动态规划,求指标函数最小值的逆序算法递归计算程序3.%输入参数4.% x各阶段状态变量的可能取值，第k列代表第k阶段状态变量可能取值5.% DecisFun(k,x)决策函数，由阶段k的状态变量x求出相应的允许决策变量6.% ObjFun(k,x,u)阶段指标函数7.% TransFun(k,x,u)状态转移函数，其中x是阶段k的某状态变量，u是相应的决策变量8.%输出参数9.% p_opt动态规划的规划过程，p_opt=[阶段序号，状态变量，决策变量，指标函数]10.% fval总目标函数值，是一个列向量，第i元素代表第一个状态变量取第i个可能值时的总目标复制代码该函数的下载参见【数学建模工具箱】/thread-221-1-2.html四、动态规划实例分析某公司拟将某种设备5台分配给甲、乙、丙3个工厂，各工厂利润与设备数量之间的关系如下表所示，问这5台设备如何分配使3个工厂的总利润为最大?应用动态规划方法分析如下1.阶段k将问题按工厂分为3个阶段，k=1,2,32.状态sk给第k个工厂分配前拥有的设备台数，显然s1=53.决策uk分配给第k个工厂的设备台数，显然分给第一个工厂可以0到s1台之间分给第二个工厂也可以0到s2台之间分给第三个工厂的为s2台1.function u=decisfun(k,s,u)2.if k==33.u=s;4.else5.u=0:s;6.end复制代码4.状态转移Tk前后两个状态之间的关系如下s2=s1-u1，s3=s2-u21.function s_next=transfun(k,s,u)2.s_next=s-u复制代码5.阶段指标Vk第k阶段的指标函数，表示配给第k个工厂uk台设备所获得的利益，显然Vk=w(uk,k)1.function V=subobjfun(k,s,u)2.w=[0 0 03. 3 5 44.7 10 65.9 11 116.12 11 127.13 11 12];8.w=-w;%由于函数只能求最小值，现在求最大值，故取符号9.%第k阶段，决策变量为u时，对应的目标值10.V=([0 1 2 3 4 5]==u)*w(:,k);%或者直接使用V=w(u,k)复制代码6.各阶段状态变量可能取值由已知，我们容易知道s1={5}s2={0,1,2,3,4,5}s3={0,1,2,3,4,5}故1.2.s=nan*ones(6,3);%没有取值的地方使用nan代替3.s(1,1)=5;4.s(:,2)=[0 1 2 3 4 5]’;5.s(:,3)=[0 1 2 3 4 5]’;复制代码根据上面的分析，我们编写程序如下1.function matlabsky2.%动态规划函数求解问题演示实例3.%by dynamic4.%see also 5.%2008.12.236.%7.%计算各状态变量可能取值，第k列代表第k个状态变量的可能取值，没有的使用NaN代替8.s=nan*ones(6,3);9.s(1,1)=5;10.s(:,2)=[0 1 2 3 4 5]';11.s(:,3)=[0 1 2 3 4 5]';12.%直接调用dynprg函数13.[p_opt,fval]=dynprog(s,@DecisFun,@ObjFun,@TransFun)14.15.function u=DecisFun(k,s,u)16.%决策函数17.if k==318. u=s;19.else20. u=0:s;21.end22.function s_next=TransFun(k,s,u)23.%状态转移函数24.s_next=s-u;25.26.function V=ObjFun(k,s,u)27.%阶段目标函数28.w=[0 0 029.3 5 430.7 10 631.9 11 1132.12 11 1233.13 11 12];34.w=-w;%由于函数只能求最小值，现在求最大值，故取符号35.%第k阶段，决策变量为u时，对应的目标值36.V=([0 1 2 3 4 5]==u)*w(:,k);%或者直接使用V=w(u,k)复制代码运行结果如下1.p_opt =2.3. 1 5 2 -74. 2 3 2 -105. 3 1 1 -46.7.8.fval =9.10. -2111. NaN12. NaN13. NaN14. NaN15. NaN复制代码这个运行结果解读如下，我们要一行一行的解读对于p_opt：第一阶段，状态变量为5，决策变量为2，阶段指标为-7，也就是说第一阶段时有5台设备，分配给第一个工厂2台，该工厂的利益为-7第二阶段，状态变量为3，决策变量为2，阶段指标为-10，也就是说第二阶段时有3台设备，分配给第二个工厂2台，该工厂的利益为-10第三阶段，状态变量为1，决策变量为1，阶段指标为-4，也就是说第三阶段时有1台设备，分配给第三个工厂1台，该工厂的利益为-4对于fval：总指标值为-21=-7-10-4，也就是说按上面的最优策略，将最大获利21后面的那些NaN是由于第一个状态变量的可能取值我们只是输入了5，而其他可能值都是用NaN的，故没有结果当然第一个状态变量的可能取多个值，dynprog可以求解第一个状态变量多取值的情况。
可是根据该题实际我们只有一个5。
下面我们试试，假如s={4 5}的运行结果1.%计算各状态变量可能取值，第k列代表第k个状态变量的可能取值，没有的使用NaN代替2.s=nan*ones(6,3);3.s(1:2,1)=[4,5];4.s(:,2)=[0 1 2 3 4 5]';5.s(:,3)=[0 1 2 3 4 5]';6.%直接调用dynprg函数7.[p_opt,fval]=dynprog(s,@DecisFun,@ObjFun,@TransFun)复制代码运行结果如下，至于结果的解读，大家可以根据dynprog函数的说明试试看看，能否明白如果没法理解，可以与我一起探讨下matlabsky@1.p_opt =2.3. 1 4 2 -74. 2 2 2 -105. 3 0 0 06. 1 5 2 -77. 2 3 2 -108. 3 1 1 -49.10.11.fval =12.13. -1714. -2115. NaN16. NaN17. NaN18. NaN书中43页的关于动态规划的编程题注：clear;clc;qita = 0;m = 0;n = 0;dE = 400*ones(3,1);dD = 400*ones(3,1);dC = 400*ones(4,1);dB = 400*ones(2,1);dA = 400*ones(1,1);F = [4;3];E = [3 5;5 2;6 6];D = [2 2 100;100 1 2;100 3 3];C = [6 8 100;3 5 100;100 3 3;100 8 4];B = [1 3 6 100;100 8 7 6];A = [5;3];% dF1 = 4;dF2 = 3; E E 到G 的最短路线disp('F1 ---> G');disp('F2 ---> G');for i = 1:length(E)for j = 1:length(F)qita = F(j) + E(i,j);if qita < dE(i)dE(i) = qita;m = i;n = j;endenddisp(['E',num2str(m),'---> F',num2str(n)]); end% Dfor i = 1:length(D)for j = 1:length(E)qita = dE(j) + D(i,j);if qita < dD(i)dD(i) = qita;m = i ;n = j;endenddisp(['D',num2str(m),'---> E',num2str(n)]); end% Cfor i = 1:length(C)for j = 1:length(D)qita = dD(j) + C(i,j);if qita < dC(i)dC(i) = qita;m = i;n = j;endenddisp(['C',num2str(m),'---> D',num2str(n)]); end% Bfor i = 1:2for j = 1:length(C)% 4qita = dC(j) + B(i,j);if qita < dB(i)dB(i) = qita;m = i;n = j;endenddisp(['B',num2str(m),'---> C',num2str(n)]);end% Afor j = 1:2qita = dB(j) + A(j);if qita < dA(1)dA(1) = qita;n = j;endenddisp(['A ---> B',num2str(n)]);disp(['A ---> G 的最短路线为：',num2str(dA)]);有用网站：/view/eb614a68a45177232f60a249.html / 全国大学生数模竞赛官方网站▲ 中国数学建模网站▲ 中科大数模网站▲/mmb/index.php 浙江大学数模网▲/hmcm 哈工大数模网站▲/~fangq/wiki/?MathTools_FAQ 数学工具FAQ▲ftp:///Public/Document/Science/matlab/▲ftp://202.198.71.195▲ftp://166.111.8.229▲ftp://166.111.30.174/forum.php（MATLAB中文论坛）。
