鐘金宏 黃 玲
1.合肥工業(yè)大學,合肥,230009
2.過程優(yōu)化與智能決策教育部重點實驗室,合肥,230009
本文研究一個實際奶制品企業(yè)的經(jīng)濟批量問題。奶制品企業(yè)為迅速進入某地市場,在當?shù)卦O立加工分廠,分廠鮮奶主要來自總部的自有規(guī)?;B(yǎng)殖場,每天由奶罐車循環(huán)運奶,分廠獨立經(jīng)營需付費給總部;總部的鮮奶產(chǎn)量和奶罐車數(shù)量有限,而客戶的需求是變動的,因此分廠有時需要從當?shù)夭少忰r奶;但為了保障品質(zhì),將優(yōu)先考慮自產(chǎn)鮮奶。根據(jù)客戶訂單或市場預測,分廠可以計算出自己的鮮奶需求;但由于奶制品銷售受天氣和節(jié)假日等影響較大,在春節(jié)等銷售旺季,經(jīng)常出現(xiàn)斷貨和延期交貨的現(xiàn)象。該問題可建模為動態(tài)離散經(jīng)濟批量模型,從總部運奶可看做是生產(chǎn),從本地采購鮮奶可看作是外包,分廠的鮮奶需求可綜合考慮生產(chǎn)、外包、庫存和延期交貨等策略來滿足。
Wagner等[1]最先研究了該問題,但未考慮外包和生產(chǎn)能力限制。生產(chǎn)能力受限的經(jīng)濟批量問題是 NP難題[2-3],僅在一些特殊情況下是多項式可解的[3-8]。外包是企業(yè)間協(xié)作的常見形式,也是企業(yè)應對復雜多變需求的主要手段。目前,考慮外包(或失去銷售)的經(jīng)濟批量問題研究很少,文獻[9-13]中關于外包、轉包和失去銷售方面的經(jīng)濟批量研究可歸為此類研究,因為它們從建模的角度來說是一樣的?,F(xiàn)有研究可分為三類:無能力限制情況[9]、庫存有界情況[10]和生產(chǎn)能力受限情況[11-13]。Sandbothe等[11]分別研究了恒定和時變生產(chǎn)能力的情況,但成本參數(shù)是常量,且未考慮失去銷售量限制。Sandbothe等[12]又研究了有時變庫存和生產(chǎn)能力限制的同樣問題,但模型參數(shù)仍為常量。Atamtürk等[13]考慮了恒定生產(chǎn)能力約束情況,但未考慮外包量的限制,分別研究了一般凹成本和非投機性成本函數(shù)情況。文獻[11-13]均沒有考慮延期交貨策略。
根據(jù)問題實際情況,本文同時考慮時變成本參數(shù)、時變生產(chǎn)能力約束、外包量受限于當前周期需求和允許延期交貨情況,提出了一個新的基于群體可行狀況和個體約束違背程度的自適應罰函數(shù),開發(fā)了相應的遺傳算法,并通過大量仿真實驗驗證了所提自適應罰函數(shù)的有效性,比較了自適應罰函數(shù)相對其他罰函數(shù)的優(yōu)勢。
通過推算,分廠可獲知在未來T周期規(guī)劃時段上的需求dt,dt>0,t=1,2,…,T。令xt和Lt分別表示周期t的生產(chǎn)量和外包量;It表示周期結束時的庫存;Ct表示周期t的生產(chǎn)能力;pt(xt)、λt(Lt)、ht(It)分別表示生產(chǎn)、外包和庫存/延期交貨成本函數(shù),它們?yōu)橐话惆己瘮?shù),則引言中描述的問題(P1)可表示為
目標函數(shù)(1)用來最小化規(guī)劃時段內(nèi)的總體生產(chǎn)、外包、庫存及延期成本;狀態(tài)等式(2)為物料平衡公式;不等式(3)表示每周期生產(chǎn)量是非負的,且不超過當前周期生產(chǎn)能力;約束(4)要求每周期外包量非負,且不超過當周期需求;等式(5)表示規(guī)劃期開始和結束時庫存為0。顯然,該問題為NP難題。
遺傳算法是求解優(yōu)化問題的有力工具,但其自身有較多參數(shù)需要確定,合理設定這些參數(shù)需根據(jù)問題具體處理,主觀性較強。根據(jù)當前和前面的搜索情況,自適應調(diào)整算法參數(shù)是一個很好的選擇,在文獻[14-19]中都有報道。其中,報道最多的是交叉及變異概率自適應和罰函數(shù)參數(shù)自適應。本文設計了一種新的能根據(jù)搜索過程反饋信息調(diào)整罰函數(shù)系數(shù)的自適應遺傳算法,圖1為算法流程圖。
圖1 自適應遺傳算法流程圖
遺傳算法本質(zhì)上只能求解無約束規(guī)劃問題,對約束規(guī)劃問題,常用處理方法有:罰函數(shù)方法;設計不可行解修正算法;設計特殊的編碼和遺傳算子。問題(P1)是約束規(guī)劃問題,存在兩類約束:每周期存在生產(chǎn)和外包能力約束;要求規(guī)劃期開始和結束時庫存為零。
在二進制編碼的遺傳算法中,可行解是針對顯型空間來說的,在基因型層次,并沒有可行和不可行說法,因此,可通過設計特殊的解碼算子,直接將基因空間的基因型映射為顯型空間的可行顯型,從而解決不可行問題。在處理第一類約束時使用了該思想,有效地降低了需處理約束的個數(shù),減小了算法在獲取可行解方面的代價。對n位二進制編碼的基因段,假定它對應的十進制數(shù)為value,其在十進制顯型空間的可行范圍為[fmin,fmax],則可采用下式將它映射入可行區(qū)域:
式中,fvalue為映射的可行顯型值。
問題假定規(guī)劃期開始時庫存為零,如不為零則可通過修正初始周期約束來轉換,成為開始庫存為零的等價問題。規(guī)劃期結束時庫存為零,這是所研究規(guī)劃問題的要求,但遺傳算法是一種有導向的隨機搜索過程,得到的解不一定能保證最終庫存為零,因此,它是算法需要處理的約束。但如因管理需要,在規(guī)劃期結束要保有庫存,這部分庫存在研究上可處理成最后周期需求。這里,規(guī)劃期結束庫存不為零是因算法本身造成的,而非管理策略要求的。
靜態(tài)罰函數(shù)法易導致算法早熟或收斂緩慢,動態(tài)罰函數(shù)法未根據(jù)算法的實際運算情況動態(tài)調(diào)整懲罰系數(shù),自適應懲罰[14-16]利用迭代過程的反饋信息動態(tài)調(diào)整懲罰系數(shù),避免了上述問題,在一些文獻的報道中均取得了較好效果。文獻[14-15]中的自適應懲罰函數(shù)方法均不同程度地需要人工確定附加參數(shù),影響了算法的可擴展性。文獻[16]設計了一個基于當前群體中不可行的狀況來改變懲罰程度的方案,避免了確定附加參數(shù)。該方法有利于避免算法早熟,但由于它過于強調(diào)不可行個體參與遺傳操作,經(jīng)常無法收斂至可行解,或搜索過程過于緩慢。特別是當群體中可行個體很少或沒有可行個體時,由于懲罰微不足道,會使搜索發(fā)散。本文測試了該方法的可用性,在12周期、24周期、36周期和48周期規(guī)劃問題上,對選擇、交叉和定標算子的120種組合分別進行實驗,每個組合運行10次,沒有一次收斂至可行解。因此,無法使用該懲罰方案解決本文問題。
為此,本文基于當前群體中可行的狀況來確定懲罰系數(shù)。在當前群體中,可行個體多時,說明搜索是收斂的,應降低懲罰程度,鼓勵不可行個體參與遺傳操作,擴大搜索區(qū)域,避免早熟收斂;反之,可行個體少時,說明搜索是發(fā)散的,應加大對可行個體懲罰,促其收斂。因此,本懲罰方法利用了群體特定信息(可行和不可行區(qū)域的比例)和染色體特定信息(不可行解距離可行區(qū)域的距離)。引入自適應懲罰后的目標函數(shù)為
其中,pop_size為群體尺寸;feasibles為當前群體中可行解個數(shù),隨GA迭代動態(tài)變化,取值范圍為0~pop_size,feasibles為0時懲罰系數(shù)取一大數(shù),最低懲罰系數(shù)是1。
下面通過仿真實驗,測試本文所提自適應罰函數(shù)的有效性。問題實例隨機生成,采用文獻[8]中的季節(jié)性需求模式,其他參數(shù)由定義在某區(qū)間上的均勻分布隨機生成(具體設定區(qū)間見表1),其中生產(chǎn)能力和需求取整數(shù),其他參數(shù)為實數(shù)。實例參數(shù)設定區(qū)間見表1。
表1 問題參數(shù)設定區(qū)間
在12周期、24周期、36周期及48周期規(guī)劃問題上,針對不同算子組合進行靜態(tài)懲罰(SPF)[17]、動 態(tài) 懲 罰 1(DPF1)[18]、動 態(tài) 懲 罰 2(DPF2)[19]、自適應懲罰1(APF1)[15]和本文所提自適應懲罰方案2(APF2)的對比實驗,從而驗證APF2的有效性,同時可確定最好的算子組合。實驗中也比較了退火懲罰[14],但其不能得到可行解,故沒有列出。
令={單點,兩點,均勻,奇偶}表示交叉算子集合,c′i(i=1,2,3,4)表示第i個交叉算子。類似地,R′={賭輪,聯(lián)賽,排序,隨機均勻采樣,確定性采樣,剩余隨機抽樣}表示選擇算子集合,r′j(j=1,2,…,6)表示第j個選擇算子,?= {無定標,線性,冪函數(shù),σ截斷,適應度值共享}表示定標方法的集合,?k(k=1,2,…,5)為第k個定標方法。c′ir′j?k表 示一個算 子 組 合,共 有120種 組合,按從左向右、深度優(yōu)先原則編號。與5種懲罰方法結合,共有600種組合,每個組合運行10次,共運行6000次。實驗中,最大迭代代數(shù)為800,群體尺寸為50,交叉概率為0.9,變異概率為0.01。分別測試了12周期、24周期、36周期及48周期規(guī)劃問題。為全面反映實驗情況,實驗結果以圖形形式提供,見圖2。圖中,橫坐標為算子組合標號(1~120);縱坐標為10次運行的最終解的成本平均值;“○”代表該算子組合的10次運行最終解均不可行的情況。
由圖2可得:①在5種罰函數(shù)中,本文提出的自適應懲罰方法(APF2)效果最好,在4個規(guī)劃問題上對所有算子組合都獲得了最低平均成本。②總體上,采用聯(lián)賽和排序選擇的組合效果較好,這兩類組合在4個規(guī)劃問題上都優(yōu)于其他組合,且對不同懲罰方法該結論基本上都成立;具體來說,對于12周期和24周期規(guī)劃問題,聯(lián)賽選擇的最終解均值要低于排序選擇;對于36周期規(guī)劃問題,排序選擇的最終解均值略好于聯(lián)賽選擇或基本相當;對于48周期問題,排序選擇的最終解均值低于聯(lián)賽選擇。
圖2 算子和懲罰方法的組合實驗結果
表2為表現(xiàn)較好的算子組合c′2r′2?4和c′1r′3?4的運算結果展示,從表2中可看出,本文提出的自適應懲罰方法效果最好,10次運行的最終解均值小于其他懲罰方案。
表2 不同懲罰方法在4個規(guī)劃問題上的表現(xiàn)
交叉和變異概率對GA的性能和行為有重要影響,常見的算子概率處理方式有:①選擇算子概率的指導規(guī)則,它未必適合所要解決的問題;②自主自適應調(diào)整方案,對算子概率進行編碼,進行同步優(yōu)化,其計算代價較大;③自適應動態(tài)調(diào)節(jié)研究,自適應動態(tài)調(diào)節(jié)在無約束優(yōu)化問題中有很多報道,在約束優(yōu)化問題上有針對可行群體的概率調(diào)整報道,但保證群體可行需要采用特殊的約束處理方法,這比較困難。本文提出的懲罰函數(shù)處理約束方法不能保證每代群體中的所有個體都可行,因此,本文采用實驗嘗試法選擇適合的算子概率,同時也觀察算法對算子概率的敏感性。
組合實驗條件為:交叉概率變動范圍為0.5~1.0,步長為0.05;變異概率變動范圍為0~0.475,步長為0.025,每個組合運行10次。實驗中,采用兩點交叉、排序選擇和適應度共享定標。圖3所示為在12周期問題上的實驗結果。圖中,橫坐標為交叉概率Pc和變異概率Pm,在間隔[b,b+0.05)內(nèi),Pc=b,Pm從 0 開 始 以 步 長0.025增至0.475;縱坐標分別為10次運行可行最終解的次數(shù)、收斂代數(shù)和成本的均值;“○”代表該組合的10次運行最終解均不可行的情況。從圖3可看出,交叉概率在0.5~1.0之間變化,算法對其不敏感;變異概率應大于零且小于0.025,在該范圍外會出現(xiàn)問題解不可行情況。后續(xù)實驗中取交叉概率為0.9,變異概率為0.01。
圖3 概率組合對最終解可行的影響
對12周期、24周期、36周期和48周期規(guī)劃問題和5種懲罰方法,算法分別運行50次。根據(jù)前面實驗,采用兩點交叉、排序選擇和適應度共享定標,交叉概率為0.9,變異概率為0.01。由于懲罰方法APF1在36周期和48周期問題上出現(xiàn)了不可行最終解,對這2個問題僅給出4種懲罰方法的運行結果。
數(shù)據(jù)表適合少量數(shù)據(jù)的比較,但對于大量數(shù)據(jù)比較,則顯得雜亂和不直觀。由于實驗中數(shù)據(jù)比較多,故采用雷達圖顯示實驗結果,見圖4。圖中,極軸標號1~50代表算法50次運行的序號;每次運行的最終解對應成本作為極徑,標注在相應的極軸上,構成一個極坐標點,將50次運行的極坐標點連起來形成一個不規(guī)則閉環(huán),描述了某種懲罰方法50次運行中最終解對應成本的變化。不同懲罰方法對應顏色深淺不同的不規(guī)則閉環(huán)。位于最內(nèi)部的不規(guī)則閉環(huán)不會被其他閉環(huán)覆蓋,表示采用該懲罰方法,最好染色體的成本最低。從圖4可看出,在4個不同周期的規(guī)劃問題上,本文提出的懲罰方法APF2對應的不規(guī)則閉環(huán)位于羅盤圖的最內(nèi)部,性能最好。
運用本文設計的自適應懲罰遺傳算法,可得到奶制品加工廠的鮮奶采購計劃,限于篇幅,這里以12周期規(guī)劃問題為例給出優(yōu)化結果,見表3。
圖4 不同懲罰方法的50次運行結果的雷達圖比較
表3 12周期規(guī)劃問題的運行結果(總成本為11862.5)
本文研究了生產(chǎn)能力受限的動態(tài)經(jīng)濟批量問題,該問題抽象于一個實際奶制品加工廠的鮮奶采購問題,涉及企業(yè)如何合理運用內(nèi)購(生產(chǎn))、外購(外包)和延期交貨等策略來最優(yōu)地滿足客戶需求。提出了一個新的自適應罰函數(shù),設計了相應的遺傳算法。與其他自適應罰函數(shù)相比,本懲罰方案不需人工設定附加參數(shù),且綜合考慮了當前群體可行狀況和個體距離可行的距離。設計了大量仿真實驗,驗證了本懲罰方案的有效性和優(yōu)勢。
[1]Wagner H M,Whitin T M.Dynamic Version of the Economic Lot-size Model[J].Management Science,1958,5(1):89-96.
[2]Florian M,Lenstra J,Kan A R.Deterministic Production Planning:Algorithms and Complexity[J].Management Science,1980,26(7):669-679.
[3]Bitran G,Yanasse H.Computational Complexity of the Capacitated Lot Size Problem[J].Management Science,1982,28(10):1174-1186.
[4]Florian M,Klein M.Deterministic Production Planning with Concave Costs and Capacity Constraints[J].Management Science,1971,18(1):12-20.
[5]Janannathan R,Rao M.A Class of Deterministic Production Planning Problems[J].Management Science,1973,19(11):1295-1300.
[6]Hoesel C P,Wagelmans A P.An O(T3)Algorithm for the Economic Lot-sizing Problem with Constant Capacities[J].Management Science,1996,42(1):142-150.
[7]Chung C,Lin C.An O(T2)Algorithm for the NI/G/NI/ND Capacitated Lot Size Problem[J].Management Science,1988,34(2):420-426.
[8]Heuvel W,Wagelmans A P.An Efficient Dynamic Programming Algorithm for a Special Case of the Capacitated Lot-sizing Problem[J].Computers &Operations Research,2006,33(12):3583-3599.
[9]Aksen D,Altinkemer K,Chand S.The Single-item Lot-sizing Problem with Immediate Lost Sales[J].Eur.J.Oper.Res.,2003,147(3):558-566.
[10]Chu F,Chu C.Single-item Dynamic Lot Sizing Models with Bounded Inventory and Outsourcing[J].IEEE Transactions on Systems,Man and Cybernetics,Part A:Systems and Humans,2008,30(1):70-77.
[11]Sandbothe R A,Thompson G L.A Forward Algorithm for the Capacitated Lot Size Model with Stockouts[J].Operation Research,1990,38(3):474-486.
[12]Sandbothe R A,Thompson G L.Decision Horizons for the Capacitated Lot Size Model with Inventory Bounds and Stockouts[J].Computers & Operations Research,1993,20(5):455-465.
[13]Atamtürk A,Hochbaum D S.Capacity Acquisition,Subcontracting,and Lot Sizing[J].Management Science,2001,47(8):1081-1100.
[14]Coello C A.Theoretical and Numerical Constraint-Handling Techniques Used with Evolutionary Algorithms:a Survey of the State of the Art[J].Computer Methods in Applied Mechanics and Engineering,2002,191(11/12):1245-1287.
[15]Hadj-Alouane A B,Bean J C.A Genetic Algorithm for the Multiple-Choice Integer Program[J].Operations Research,1997,45(1):92-101.
[16]Dadios E P,Ashraf J.Genetic Algorithm with Adaptive and Dynamic Penalty Functions for the Selection of Cleaner Production Measures:A Constrained Optimization Problem[J].Clean Tech.Environ Policy,2006,8(2):85-95.
[17]Homaifar A,Qi C X,Lai S H.Constrained Optimization Via Genetic Algorithms[J].Simulation 1994,62(4):242-254.
[18]Joines J,Houck C R.On the Use of Non-stationary Penalty Functions to Solve Nonlinear Constrained Optimization Problems with Gas[C]//Proceedings of the First IEEE International Conference on Evolutionary Computation.Orlando,1994:579-584.
[19]Kazarlis S,Petridis V.Varying Fitness Functions in Genetic Algorithms:Studying the Rate of Increase in the Dynamic Penalty Terms[C]//Proceedings of the 5th Int.Conf.on Parallel Problem Solving from Nature.Berlin,1998:211-220.