1、第14章 動態規劃本章通過將實際問題當作多階段決策過程來建立數學模型，學習和掌握動態規劃的相關知識及其求解的逆序算法，并結合計算機編程解決實際問題。14.1 引例：生產計劃的制定某工廠與用戶訂立合同，在四個月內出售一定數量的某種產品，產量限制為10的倍數，工廠每月最多生產100件，產品可以存儲，存儲費用為每臺2百元，每個月的需求量及每件產品的生產成本見表14.1。表14.1 生產成本和需要量月份每件生產成本（元）需要量（件）170602727038012047660現在分別在（1）一月初沒有存貨可用和（2）一月初有20件存貨可用這兩種情況下確定每月的出產量，要求既能滿足每月的合同需求量，又使生

2、產成本和存儲費用達到最小。靜態的看，本引例是一個整數規劃問題，但這里我們可以把這個問題的解決動態地視為各月（一般稱為階段）先后作出決策（這里指生產量）的過程多階段的決策過程，而在每個月作決策時，不能僅考慮本月的費用（一般稱為階段指標）因為本月的決策會對以后各月的決策產生影響，因此應考慮從本月直到第四月末的總費用（總指標），而每月的決策依賴于各月初倉庫的存貨量（一般稱為始端）和以前各月如何造成這貨存量的情況無關（稱為無后效性）。在例2中我們將計算一月初無存貨可用時的最優決策見表2表14.2 最優生產計劃1月2月3月4月月初存儲量（件）040700產量（件）1001005060則第四月的決策即為月

3、初倉儲數為0時的最優決策，第三、四月的決策即為第三月初倉儲數為70時的最優決策，以及第二、三、四月的決策即為第二月初倉儲數為40時的最優決策。因為若不然，如對應于第三月初倉儲數為70時，第三、四月的最優決策是分別生產80件和30件（即這樣的費用比分別生產50件和60件更?。?，則我們保留第一、二月的生產數，而把第三第四月分別改為80件和30件，這個方案顯然優于原來的方案，這和原來的方案是最優相矛盾。這個性質可以簡述為：最優決策的任何截斷仍是最優的（最優性原理）。把這一最優化問題視為符合最優性原理、無后效性的多階段決策過程并進行求解的方法稱為動態規劃方法。14.2 動態規劃的基本理論復習14.2.

4、1基本思想與逆序解法的直觀回顧前面我們已簡單介紹了動態規劃。為了更便于 了解動態規劃的基本思想、描述方式和逆序解法，我們來看一個確定網絡最短路徑問題的例子。例1：（最短路徑的確定）以下是一個賦權圖，兩頂點連線上的數字表示距離，確定一條從始點到終點鋪設管道并使總距離最短的路線。 V4 1 6 2 V11 3 V2 3 3 V8 2 5 V14 4 5 6 V5 5 8 1 V12 2 5 3 V16V1 3 V9 2 6 6 V15 3 8 7 V6 3 V13 V3 6 8 3 V10 3 V7 4 圖14.1直觀上我們有這樣一個重要常識：如果由起點經過點和點而到達終點是一條最短路線，則由點出

5、發經過點到達終點的這條路線，對于從點出發 到達終點的所有可能選擇的不同路線來說，必定也是最短路線。例如在最短路線問題中，若找到了是由點到點最短路線，則應該是由點出發到達終點的所有可能選擇的不同路線中最短路線。這一特征即為上面所提到的最優性原理：最優性決策的任何截斷仍是最優的，這是動態規劃的基本原理。根據最優性原理，尋找最短路線可從最后一段開始，用由后向前逐步遞推的方法，求出各點到點的最短路線，最后求得由點到點的最短路線，所以動態規劃的逆序求解方法是從終點逐段向始點方向尋找最短路線的一種方法。下面逐段完成計算。例1當然可以用圖與網絡實驗中介紹的Dijkstra算法來求解，但這里我們將把該問題看成

6、6個階段的決策過程，并逆序逐段求解，從而較直觀地揭示動態規劃的基本思想。時，以表示由到的最短距離，以表示由到的最短距離，則，。時，出發點有三個、和。若從出發有和兩個選擇，以表示由到的最短距離，表示由到的距離，表示由到的距離，表示相應的選擇或決策，則：可見，其最短路徑為。同理，從出發也有和兩個選擇，、和意義與上面相似，則：可見，其最短路徑為。從出發，同樣有：可見，其最短路徑為。時，有三個出發點、和，同樣計算如下：可見，其最短路徑為?？梢?，其最短路徑為。可見，其最短路徑為。時，同樣計算有：，；，；，時，同樣計算有：，；，時，只有一個出發點，則：可見，所以本題的最短路徑為。14.2.2 動態規劃的基

7、本概念及其數學描述（1）階段 整個問題的解決可分為若干個相互聯系的階段依次進行。通常按時間或空間劃分階段，描述階段的變量稱為階段變量，計為。（2）狀態 狀態表示每個階段開始所處的自然狀況或客觀條件，它描述了研究問題過程的狀況。個階段的狀態通常用狀態變量描述。常用表示第階段的狀態變量。個階段的決策過程有個狀態。用規劃方法解決多階段決策問題時，要求整個過程具有無后效性，即：如果某階段的狀態給定，則此階段以后過程的發展不受以前狀態的影響，未來狀態只依賴于當前狀態。（3）決策 某一階段的狀態確定后，可以作出各種選擇從而演變到下一階段某一狀態，這種選擇手段稱為決策。描述決策的變量稱為決策變量。決策變量限

8、制的取值范圍稱為允許決策集合。用表示第階段處于狀態時的決策變量，它是的函數，用表示的允許決策的集合。比如在例1中，。（4）策略 一個由每個階段的決策按順序排列組組成的集合稱為策略，用表示。即。由第階段的狀態開始到終止狀態的后部子過程的策略記為，即。在實際問題中，可供選擇的策略有一定范圍，此范圍稱為允許策略集合。允許策略集合中達到最優效果的策略稱為最優策略。（5）狀態轉移方程 如果第個階段狀態變量為，作出的決策，那么第階段的狀態變量也被完全確定。用狀態轉移方程表示這種演變規律，寫作。（6）指標函數和最優值函數 指標函數是系統執行某一策略所產生結果的數量表示，是用來衡量策略優劣的數量指標，它定義在

9、全過程和所有后部子過程上，分別用和表示。即和。過程在某個階段的階段指標函數（或階段效益）是衡量該階段決策優劣的數量指標，它取決于狀態和決策，用表示。如例1中兩頂點間的距離是階段指標函數，而指標函數則是直到終點的距離的和。根據不同的實際問題，效益可以是利潤、距離、產量或資源等。指標函數往往是各階段效益的某種形式。指標函數的最優值稱為最優函數。（7）最優策略和最優軌線 使指標函數達到最優值的策略是從階段開始的后部子過程的最優策略，記為。是全過程的最優策略，簡稱為最優策略。從初始狀態出發，過程按照和狀態轉移方程演變所經歷的狀態序列稱為最優軌線。14.3 動態規劃逆序算法的MATLAB程序14.3.1

10、逆序算法的基本方程由例1的求解過程可以看出下面的方程在動態規劃逆序求解中起著本質的作用稱此為動態規劃逆序求解的基本方程?？梢园呀討B規劃模型歸納成以下幾個步驟（1）將問題恰當地劃分為若干個階段；（2）正確選擇狀態變量，使它既能描述過程的演變，又滿足無后效性；（3）規定決策變量，確定每個階段允許決策集合；（4）寫出狀態轉移方程；（5）確定個階段各種決策的階段指標，列出計算各階段最優后部策略指標的基本方程。14.3.2動態規劃逆序算法的MATLAB程序Dynprog.m%M函數dynprog.mfunctionp_opt,fval=dynprog(x,DecisFun,ObjFun,TransF

11、un)% p_opt,fval =dynprog(x,DecisFun,ObjFun,TransFun)% 自由始端和終端的動態規劃，求指標函數最小值的逆序算法遞歸計算程序。x是狀態變量，一列代表一個階段狀態；% M函數DecisFun（k,x）由階段k的狀態變量x求出相應的允許決策變量；% M函數ObjFun(k,x,u)是階段指標函數，M函數TransFun(k,x,u)是狀態轉移函數，其中x 是階段k的某狀態變量，u是相應的決策變量；% 輸出p_opt由4列構成，p_opt=序號組；最優軌線組；最優策略組；% 指標函數組;fval是一個列向量，各元素分別表示p_opt各最優策略組對應始端

12、狀態x的最優函數組；k=length(x(1,:); f_opt=nan*ones(size(x); d_opt=f_opt;t_vubm=inf*ones(size(x); x_isnan=isnan(x); t_vub=inf;% 計算終端相關值tmpl=find(x_isnan(:,k);tmp2=length(tmpl);for I=1:tmp2 u=feval(DecisFun,k,x(i,k); tmp3=length(u); for j=1:tmp3 tmp=feval(ObjFun,k,x(tmpl(i),k), u(j); if tmp=t_vub, f_opt(i,k)=t

13、mp; d_opt(i,k)=u(j); t_vub=tmp;end;end;end%逆推計算各階段的遞歸調用程序for ii=k-1:-1:1 tmp10=find(x_isnan(:,ii);tmp20=length(tmp10); for i=1:tmp20 u=feval(DecisFun,ii,x(i,ii);tmp30=length(u); for j=1:tmp30 tmp00=feval(ObjFun,ii,x(tmp10(i),ii),u(j); tmp40=feval(TransFun,ii,x(tmp10(i),ii),u(j); tmp50=x(:,ii+1)-tmp4

14、0; tmp60=find(tmp50= =0); ifisempty(tmp60), tmp00=tmp00+f_opt(tmp60(1),ii+1); if tmp00=t_vubm(i,ii) f_opt(i,ii)=tmp00; d_opt(i,ii)=u(j); t_vubm(i,ii)=tmp00;end;end;end;end;end;fval=f_opt(tmp1,1);% 記錄最優決策、最優軌線和相應指標函數值p_opt=;tmpx=;tmpd=;tmpf=;tmp0=find(x_isnan(;,1);tmp01=length(tmp0);for i=1:tmp01, tm

15、pd(i)=d_opt(tmp0(i),1); tmpx(i)=x(tmp0(i),1); tmpf(i)=feval(ObjFun,1,tmpx(i),tmpd(i); p_opt(k*(I-1)+1,1,2,3,4)=1,tmpx(i),tmpd(i),tmpf(i); for ii=2:k tmpx(i)=feval(TransFun,ii-1,tmpx(i),tmpd(i); tmp1=x(:,ii)-tmpx(i); tmp2=find(tmp1= =0); ifisempty(tmp2) tmpd(i)=d_opt(tmp2(1),ii); end; tmpf(i)=feval(O

16、bjFun,ii,tmpx(i),tmpd(i); p_opt(k*(i-1)+ii,1,2,3,4)=ii,tmpx(i),tmpd(i),tmpf(i);end;end;14.4 例題下面將就動態規劃的幾個典型應用問題分別舉例計算。例2：（生產計劃制定）解：這是一個4階段動態規劃問題。如果用逆序法解題，第1階段是1月份，第4階段是4月份。記第階段開始的產品存儲數（狀態變量）；第階段的產量（決策變量）；第階段每件產品的生產成本；第階段的需求量；階段指標的函數為：；狀態轉移方程為：；基本方程為：對于本例的問題（1）：一月初無存貨，可首先分析出各月的最大存貨量和產量。如對4月初，前三個月的最大產

17、量為300件（每月最大產量為100件），實際需求量為60+70+120=250（件）。所以四月初的最大存儲量為300-250=50（件）因產量限制為10的倍數，4月初的存貨量只可能是0、10、20、30、40、50、這六種。而4月份的實際需要量為60件，因此第一階段的產量（決策）相應為60、50、40、30、20、10，進而計算，再分析3月初情況等等，如此可仿例1逐段手工計算求出最優決策。下面將調用參考程序dynprog.m進行計算。由于計算機的優勢，這里將就1月初存貨分別為0、10、20、30、40、50、60的所有可能情況進行計算。把此問題作為自由始端動態規劃考慮，此時個階段的最大存貨出現

18、在時，經分析可知，考慮到產量是10的倍數，根據上面所述的階段指標函數、狀態轉移方程和基本方程，寫出下面的3個M函數以備計算時調用，函數意義見函數的說明部分。%M函數eg13f1_2.mfunction u=DecisF_1(k,x)% 在階段k由狀態變量x的值求出相應的決策變量所有的取值c=70,72,80,76;q=10*6,7,12,6;if q(k)-xclear;x=nan*ones(14,4); %x是10的倍數，最大范圍,% 因此x=0,1,，13,所以x初始化取14行，nan表示無意義元素x(1:7,1)=10*(0:6); %按月定義x的可能取值x(1:11,2)=10*(0:

19、10); x(1:12,3)=10*(2:13);x(1:7,4)=10*(0:6);p,f=dynprog(x,eg13f1_2,eg13f2_2,eg13f3_2)p= %輸出結果1 0 100 70002 40 100 72803 70 50 41404 0 60 45601 10 100 7020 2 50 100 73003 80 40 33604 0 60 45601 20 100 70402 60 100 73203 90 30 25804 0 60 45601 30 100 70602 70 100 73403 100 20 18004 0 60 45601 40 100 70

20、80 2 80 100 73603 110 10 10204 0 60 45601 50 100 71002 90 100 73803 120 0 2404 0 60 45601 60 100 71202 100 100 74003 130 0 2604 10 50 3820f=22980222402150020760200201928018600由p 第一和第三個4行可以看出一月初無存貨和有存貨20件的最優決策，現將其列成下表，以供對比理解。表14.31月初的存貨量階段序號最優軌線(存儲量)最優決策(產量)階段指標函數值(成本)0(件)1010070002401007280370504140

21、4060456020(件)1201007040260100732039030258040604560由f的第1、3行可以看出對應四個月的總成本分別為22980和21500。例3：調用dymprog.m計算例1中網絡的最短路徑。解：首先編寫3個M函數。以供計算時調用。由表1可知狀態變量有7個，并可寫出下面的函數。%M函數eg13f1_3.mfunction u=eg13f1_3(k,x)% 在階段k由狀態變量x的值求出其相應的決策變量所有的取值% 這里求出了所有狀態x對應的決策變量值uif x= =1,u=2;3;elseif x= =2,u=4;5;6;elseif x= =3,u=5;6;7

22、;elseif(x= =4)|(x= =5),u=8;9;elseif(x= =6)|(x= =7),u=9;10; elseif x= =8,u=11,12;elseif(x= =9)|(x= =10),u=12;13;elseif(x= =11)|(x= =12)|(x= =13),u=14,15;elseif(x= =14)|(x= =15),u=16;elseif x= =16,u=16;end %M函數eg13f2_3.mfunction v=eg13f2_3(k,x,u)% 各階段指標函數值，如x=1,u=2時v=5。tt=5;3;1;3;6;8;7;6;6;8;3;5;3;3;8

23、;4;2;2;1;2;3;3;3;3;5;5;2;6;6;4;3;tmp=x= =1&u= =2,x= =1&u= =3,x= =2&u= =4,x= =2&u= =5,x= =2&u= =6, x= =3&u= =5, x= =3&u= =6, x= =3&u= =7,x= =4&u= =8, x= =4&u= =9, x= =5&u= =8, x= =5&u= =9,x= =6&u= =9,x= =6&u= =10, x= =7&u= =9, x= =7&u= =10,x= =8&u= =11, x= =8&u= =12, x= =9&u= =12, x= =9&u= =13,x= =10

24、&u= =12, x= =10&u= =13, x= =11&u= =14, x= =11&u= =15,x= =12&u= =14, x= =12&u= =15, x= =13&u= =14, x= =13&u= =15,x= =14&u= =16, x= =15&u= =16;v=tmp*tt;%M函數eg13f3_3.mfunction y= eg13f3_3(k,x,u)y=u; %狀態轉移方程計算如下：clear; x=nan*ones(4,7); % 初始化，nan為無意義元素 % 7個狀態變量，每個最多4種取值x(1,1)=1;x(1:2,2)=2;3; %逐列定義x,其值為頂點

25、的序號x(1:4,3)=(4:7);x(1:3,4)=(8:10);x(1:3,5)=(11:13);x(1:2,6)=14;15; x(1,7)=16;p,f=dynprog(x,eg13f1_3,eg13f2_3, eg13f3_3)p= %輸出結果 1 1 2 52 2 5 33 5 8 34 8 12 25 12 15 2 6 15 16 37 16 16 0f=18 %最短路徑距離可見最短路徑按頂點序號為1 2 5 8 12 15 16。例4：某工業部門根據國家計劃的安排，擬將5臺某種高效率的設備，分配給所屬的甲、乙、丙三個工廠，各工廠若獲得這種設備之后，可以為國家提供的盈利如下面的

26、表。問：這5臺設備如何分配給各工廠才能使國家得到的盈利最大？解：將問題按工廠分為三個階段，甲、乙和丙3個工廠分別編號為1、2和3。設狀態變量表示分配給第個工廠至第個工廠的設備臺數。決策變量表示分配給第個工廠的設備臺數。則狀態轉移方程為分配給第個工廠至第個工廠的設備臺數。階段指標函數表示臺設備分配到第個工廠所獲得的盈利值。表示臺設備分配給第個工廠至第個工廠所獲得的最大盈利值。則基本方程為表14.4設備數甲乙丙000013542710639111141211125131112利用計算機計算的優勢，可根據上表將本問題視為自由始端，即初始狀態的動態規劃問題求解。由此可以寫出下面3個M函數，以供計算時調用。%M函數eg13f1_4.mfunction u=eg13f